diff --git a/doc/src/fix_plumed.html b/doc/src/fix_plumed.html new file mode 100644 index 0000000000..4342403e84 --- /dev/null +++ b/doc/src/fix_plumed.html @@ -0,0 +1,118 @@ + +
LAMMPS WWW Site - LAMMPS Documentation - LAMMPS Commands +
+ + + + + + +
+ +

fix plumed command +

+

Syntax: +

+
fix ID group-ID plumed keyword value ...
+
+ +

Examples: +

+
+  fix pl all plumed all plumed plumedfile plumed.dat outfile p.log 
+
+

Description: +

+

+This fix instructs LAMMPS to call the PLUMED library, which allows one +to perform various forms of trajectory analysis on the fly and to also use +methods such as umbrella sampling and metadynamics to enhance the sampling of +phase space. +

+

The documentation included here only describes the fix plumed command. This command +is LAMMPS specific whereas most of the functionality implemented in PLUMED will work with a +range of MD codes and also when PLUMED is used as a stand alone code. The full documentation +for PLUMED is available at this website +

+

The PLUMED library is developed at https://github.com/plumed/plumed2 +A detailed discussion of the code can be found in (PLUMED). +

+

There are some example scripts for using this package with LAMMPS in the +examples/USER/plumed directory. +

+
+ +

The command to call PLUMED above is reasonably self explanatory. Within the input file +for lammps the user is required to specify the input file for PLUMED and a file on which +to output the PLUMED log. The user must specify both of these arguments every time +PLUMED is to be used. Furthermore, the fix plumed command should appear in the LAMMPS input +file after the relevant input paramters (e.g. the timestep) have been set. +

+

The group-ID entry is ignored. LAMMPS will always pass all the atoms to PLUMED +and there can only be one instance of the plumed fix at a time. The plumed fix communicates +the minimum amount of information required and the PLUMED supports multiple, completely +independent collective variables, multiple independent biases and multiple independent forms of analysis. +There is thus really no restriction in functionality by only allowing only one plumed fix in the LAMMPS input. +

+

The plumedfile keyword allows the user to specify the name of the PLUMED input file. +Instructions as to what should be included in a plumed input file can be found in the +documentation for PLUMED. +

+

The outfile keyword allows the user to specify the name of a file on which to output +the PLUMED log. This log file normally just parots the information that is contained in the input +file. The names of the files on which the results from the various analyses that have been performed +using PLUMED will be specified by the user in the PLUMED input file. +

+

Restart, fix_modify, output, run start/stop, minimize info: +

+

+When performing a restart of a calculation that involves PLUMED you must include a RESTART command +in the PLUMED input file as detailed in the PLUMED documentation. When the restart +command is found in the PLUMED input PLUMED will append to the files that were generated in the run +that was performed previously. Furthermore, any history dependent bias potentials that were accumulated in +previous calculations will be read in when the restart command is included in the PLUMED input. +

+

The fix_modify energy option is not supported by this fix. +

+

Nothing is computed by this fix that can be accessed by any of the +output commands within LAMMPS. All the quantities +of interest can be output by commands that are native to PLUMED, however. +

+

Restrictions: +

+

This fix is part of the USER-PLUMED package. It is only enabled if +LAMMPS was built with that package. See the Making +LAMMPS section for more info. +

+

There can only be one plumed fix active at a time. Since the interface +communicates only the minimum amount of information and since the PLUMED module +itself can handle an arbitrary number of analysis and biasing methods, this is +not a limitation of functionality. +

+

Related commands: +

+

fix smd +

fix colvars +

+

Default: +

+

The default options are plumedfile = NULL and outfile = NULL +

+
+ + + +

(PLUMED) G.A. Tribello, M. Bonomi, D. Branduardi, C. Camilloni and G. Bussi, Comp. Phys. Comm 185, 604 (2014) DOI:https://doi.org/10.1016/j.cpc.2013.09.018 +

+ + diff --git a/src/.gitignore b/src/.gitignore index 92933ce5ee..b84ab65787 100644 --- a/src/.gitignore +++ b/src/.gitignore @@ -53,6 +53,11 @@ /colvarproxy_lammps_version.h /fix_colvars.cpp /fix_colvars.h +/fix_plumed.cpp +/fix_plumed.h +/Plumed.cpp +/Plumed.h +/Plumed.inc /dump_molfile.cpp /dump_molfile.h /molfile_interface.cpp diff --git a/src/USER-PLUMED/Install.sh b/src/USER-PLUMED/Install.sh new file mode 100755 index 0000000000..3e0bef56ca --- /dev/null +++ b/src/USER-PLUMED/Install.sh @@ -0,0 +1,39 @@ +# Install/unInstall package files in LAMMPS +# edit 2 Makefile.package files to include/exclude ATC info + +if (test $1 = 1) then + + if (test -e ../Makefile.package) then + sed -i -e 's|^PKG_LIB =[ \t]*|& $(PLUMED_LOAD) |' ../Makefile.package + fi + + if (test -e ../Makefile.package.settings) then + # multiline form needed for BSD sed on Macs + sed -i -e '4 i \ +include ..\/Plumed.inc +' ../Makefile.package.settings + fi + + cp fix_plumed.cpp .. + cp fix_plumed.h .. + cp Plumed.h .. + cp Plumed.cpp .. + cp Plumed.inc .. + +elif (test $1 = 0) then + + if (test -e ../Makefile.package) then + sed -i -e 's/[^ \t]* \$(PLUMED_LOAD)[^ \t]* //' ../Makefile.package + fi + + if (test -e ../Makefile.package.settings) then + sed -i -e '/^include.*Plumed\.inc.*$/d' ../Makefile.package.settings + fi + + rm -f ../fix_plumed.cpp + rm -f ../fix_plumed.h + rm -f ../Plumed.h + rm -f ../Plumed.cpp + rm -f ../Plumed.inc + +fi diff --git a/src/USER-PLUMED/Plumed.cpp b/src/USER-PLUMED/Plumed.cpp new file mode 100644 index 0000000000..f279c4a0eb --- /dev/null +++ b/src/USER-PLUMED/Plumed.cpp @@ -0,0 +1,471 @@ +/* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ + Copyright (c) 2011-2018 The plumed team + (see the PEOPLE file at the root of the distribution for a list of names) + + See http://www.plumed.org for more information. + + This file is part of plumed, version 2. + + plumed is free software: you can redistribute it and/or modify + it under the terms of the GNU Lesser General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + plumed is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public License + along with plumed. If not, see . ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */ +#include "Plumed.h" + +#ifdef __PLUMED_HAS_DLOPEN +#include +#endif +#include +#include +#include +#include +#include + +/* DECLARATION USED ONLY IN THIS FILE */ + +#ifdef __cplusplus +extern "C" { +#endif + +/** + Function pointer to plumed_create +*/ + +typedef void*(*plumed_create_pointer)(void); +/** + Function pointer to plumed_cmd +*/ +typedef void(*plumed_cmd_pointer)(void*,const char*,const void*); + +/** + Function pointer to plumed_finalize +*/ +typedef void(*plumed_finalize_pointer)(void*); + +/** + Holder for plumedmain function pointers. +*/ +typedef struct { + plumed_create_pointer create; + plumed_cmd_pointer cmd; + plumed_finalize_pointer finalize; +} plumed_plumedmain_function_holder; + +/** + Holder for plumed symbol table. +*/ +typedef struct { + int version; + plumed_plumedmain_function_holder functions; +} plumed_symbol_table_type; + +/** + Register for plumedmain function pointers +*/ +plumed_plumedmain_function_holder* plumed_kernel_register(const plumed_plumedmain_function_holder*); + +#ifdef __PLUMED_STATIC_KERNEL +/* Real interface */ +void*plumed_plumedmain_create(void); +void plumed_plumedmain_cmd(void*,const char*,const void*); +void plumed_plumedmain_finalize(void*); +#else +/* dummy interface */ +void*plumed_dummy_create(void); +void plumed_dummy_cmd(void*,const char*,const void*); +void plumed_dummy_finalize(void*); +#endif + +#ifdef __cplusplus +} +#endif + +/* END OF DECLARATION USED ONLY IN THIS FILE */ + +/* These are the dummy routines which are used when plumed is not available */ + +#ifdef __PLUMED_STATIC_KERNEL + +static int installed=1; + +#else + +static int installed=0; + +static int dummy; + +void*plumed_dummy_create(void) { + return (void*)&dummy; +} + +void plumed_dummy_cmd(void*p,const char*key,const void*val) { + (void) p; /* avoid warning on unused parameter */ + (void) key; /* avoid warning on unused parameter */ + (void) val; /* avoid warning on unused parameter */ + fprintf(stderr,"+++ ERROR: you are trying to use plumed, but it is not available +++\n"); + fprintf(stderr,"+++ Check your PLUMED_KERNEL environment variable +++\n"); + exit(1); +} + +void plumed_dummy_finalize(void*p) { + (void) p; /* avoid warning on unused parameter */ +} + +#endif + +plumed_plumedmain_function_holder* plumed_kernel_register(const plumed_plumedmain_function_holder* f) { + /* + Argument f is present for historical reasons but ignored in PLUMED>=2.5. + */ + if(f) { + if(getenv("PLUMED_LOAD_DEBUG")) { + fprintf(stderr,"+++ Ignoring registration at %p (%p,%p,%p) +++\n",(void*)f,(void*)f->create,(void*)f->cmd,(void*)f->finalize); + } + } +#ifdef __PLUMED_STATIC_KERNEL + /* + When __PLUMED_STATIC_KERNEL is defined, the function holder is initialized + to statically bound plumed_plumedmain_create, plumed_plumedmain_cmd, plumed_plumedmain_finalize and + cannot be changed. This saves from mis-set values for PLUMED_KERNEL. + */ + static plumed_plumedmain_function_holder g= {plumed_plumedmain_create,plumed_plumedmain_cmd,plumed_plumedmain_finalize}; +#else + /* + On the other hand, for runtime binding, we use dlsym to find the relevant functions. + Notice that as of PLUMED 2.5 self registration of the kernel is ignored, so argument f + is not used anymore. + Also notice that we should put some guard here for safe multithread calculations. + */ + static plumed_plumedmain_function_holder g= {plumed_dummy_create,plumed_dummy_cmd,plumed_dummy_finalize}; + static int first=1; +#ifdef __PLUMED_HAS_DLOPEN + const char* path; + char* pathcopy; + void* p; + char* pc; + plumed_symbol_table_type* plumed_symbol_table_ptr; + plumed_plumedmain_function_holder functions; + char* debug; + size_t strlenpath; + int dlopenmode; + /* + f==NULL is required here otherwise we would enter this block a second time + when plumed_kernel_register is called by the just loaded shared library. + */ + if(first && f==NULL) { + path=getenv("PLUMED_KERNEL"); + debug=getenv("PLUMED_LOAD_DEBUG"); +#ifdef __PLUMED_DEFAULT_KERNEL + /* + This variable allows a default path for the kernel to be hardcoded. + Can be useful for hardcoding the predefined plumed location + still allowing the user to override this choice setting PLUMED_KERNEL. + The path should be chosen at compile time adding e.g. + -D__PLUMED_DEFAULT_KERNEL=/opt/local/lib/libplumed.dylib + */ + /* This is required to add quotes */ +#define PLUMED_QUOTE_DIRECT(name) #name +#define PLUMED_QUOTE(macro) PLUMED_QUOTE_DIRECT(macro) + if(! (path && (*path) )) path=PLUMED_QUOTE(__PLUMED_DEFAULT_KERNEL); +#endif + if(path && (*path)) { + fprintf(stderr,"+++ Loading the PLUMED kernel runtime +++\n"); + fprintf(stderr,"+++ PLUMED_KERNEL=\"%s\" +++\n",path); + if(getenv("PLUMED_LOAD_NAMESPACE") && !strcmp(getenv("PLUMED_LOAD_NAMESPACE"),"LOCAL")) { + dlopenmode=RTLD_NOW|RTLD_LOCAL; + if(debug) fprintf(stderr,"+++ Loading with mode RTLD_NOW|RTLD_LOCAL +++\n"); + } else { + dlopenmode=RTLD_NOW|RTLD_GLOBAL; + if(debug) fprintf(stderr,"+++ Loading with mode RTLD_NOW|RTLD_GLOBAL +++\n"); + } + p=dlopen(path,dlopenmode); + if(!p) { + /* + Something went wrong. We try to remove "Kernel" string from the PLUMED_KERNEL variable + and load directly the shared library. Notice that this particular path is only expected + to be necessary when using PLUMED<=2.4 and the symbols in the main executable are + not visible. All the other cases (either PLUMED>=2.5 or symbols in the main executable visible) + should work correctly without entering here. + */ + fprintf(stderr,"+++ An error occurred. Message from dlopen(): %s +++\n",dlerror()); + strlenpath=strlen(path); + pathcopy=(char*) malloc(strlenpath+1); + strncpy(pathcopy,path,strlenpath+1); + pc=pathcopy+strlenpath-6; + while(pc>=pathcopy && memcmp(pc,"Kernel",6)) pc--; + if(pc>=pathcopy) { + memmove(pc, pc+6, strlen(pc)-5); + fprintf(stderr,"+++ Trying %s +++\n",pathcopy); + p=dlopen(pathcopy,dlopenmode); + if(!p) fprintf(stderr,"+++ An error occurred. Message from dlopen(): %s +++\n",dlerror()); + } + free(pathcopy); + } + if(p) { + functions.create=NULL; + functions.cmd=NULL; + functions.finalize=NULL; + /* + If the library was loaded, use dlsym to initialize pointers. + Notice that as of PLUMED 2.5 we ignore self registrations. + Pointers are searched in the form of a single pointer to a structure, which + is the standard way in PLUMED 2.5, as well as using alternative names used in + PLUMED 2.0 to 2.4 (e.g. plumedmain_create) and in some intermediate versions between + PLUMED 2.4 and 2.5 (e.g. plumed_plumedmain_create). The last chance is probably + unnecessary and might be removed at some point. + */ + plumed_symbol_table_ptr=(plumed_symbol_table_type*) dlsym(p,"plumed_symbol_table"); + if(plumed_symbol_table_ptr) functions=plumed_symbol_table_ptr->functions; + if(debug && plumed_symbol_table_ptr) { + fprintf(stderr,"+++ plumed_symbol_table version %i found at %p +++\n",plumed_symbol_table_ptr->version,(void*)plumed_symbol_table_ptr); + fprintf(stderr,"+++ plumed_function_pointers found at %p (%p,%p,%p) +++\n",(void*)&plumed_symbol_table_ptr->functions,(void*)functions.create,(void*)functions.cmd,(void*)functions.finalize); + } + + if(!functions.create) { + functions.create=(plumed_create_pointer) dlsym(p,"plumedmain_create"); + if(debug && functions.create) fprintf(stderr,"+++ %s found at %p +++\n","plumedmain_create",(void*)functions.create); + } + if(!functions.create) { + functions.create=(plumed_create_pointer) dlsym(p,"plumed_plumedmain_create"); + if(debug && functions.create) fprintf(stderr,"+++ %s found at %p +++\n","plumed_plumedmain_create",(void*)functions.create); + } + + if(!functions.cmd) { + functions.cmd=(plumed_cmd_pointer) dlsym(p,"plumedmain_cmd"); + if(debug && functions.cmd) fprintf(stderr,"+++ %s found at %p +++\n","plumedmain_cmd",(void*)functions.cmd); + } + if(!functions.cmd) { + functions.cmd=(plumed_cmd_pointer) dlsym(p,"plumed_plumedmain_cmd"); + if(debug && functions.cmd) fprintf(stderr,"+++ %s found at %p +++\n","plumed_plumedmain_cmd",(void*)functions.cmd); + } + + if(!functions.finalize) { + functions.finalize=(plumed_finalize_pointer) dlsym(p,"plumedmain_finalize"); + if(debug && functions.finalize) fprintf(stderr,"+++ %s found at %p +++\n","plumedmain_finalize",(void*)functions.finalize); + } + if(!functions.finalize) { + functions.finalize=(plumed_finalize_pointer) dlsym(p,"plumed_plumedmain_finalize"); + if(debug && functions.finalize) fprintf(stderr,"+++ %s found at %p +++\n","plumed_plumedmain_finalize",(void*)functions.finalize); + } + + if(functions.create && functions.cmd && functions.finalize) { + g=functions; + installed=1; + } else { + if(!functions.create) fprintf(stderr,"+++ pointer to (plumed_)plumedmain_create not found +++\n"); + if(!functions.cmd) fprintf(stderr,"+++ pointer to (plumed_)plumedmain_cmd not found +++\n"); + if(!functions.finalize) fprintf(stderr,"+++ pointer to (plumed_)plumedmain_finalize not found +++\n"); + } + } + } + first=0; + } +#endif +#endif + return &g; +} + +/* C wrappers: */ + +plumed plumed_create(void) { + plumed p; + plumed_plumedmain_function_holder*h=plumed_kernel_register(NULL); + assert(h); + assert(h->create); + p.p=(*(h->create))(); + assert(p.p); + return p; +} + +void plumed_cmd(plumed p,const char*key,const void*val) { + plumed_plumedmain_function_holder*h=plumed_kernel_register(NULL); + assert(p.p); + assert(h); + assert(h->cmd); + (*(h->cmd))(p.p,key,val); +} + +void plumed_finalize(plumed p) { + plumed_plumedmain_function_holder*h=plumed_kernel_register(NULL); + assert(p.p); + assert(h); + assert(h->finalize); + (*(h->finalize))(p.p); + p.p=NULL; +} + +int plumed_installed(void) { + plumed_kernel_register(NULL); + return installed; +} + +/* we declare a Plumed_g_main object here, in such a way that it is always available */ + +static plumed gmain= {NULL}; + +plumed plumed_global(void) { + return gmain; +} + +void plumed_gcreate(void) { + assert(gmain.p==NULL); + gmain=plumed_create(); +} + +void plumed_gcmd(const char*key,const void*val) { + assert(gmain.p); + plumed_cmd(gmain,key,val); +} + +void plumed_gfinalize(void) { + assert(gmain.p); + plumed_finalize(gmain); + gmain.p=NULL; +} + +int plumed_ginitialized(void) { + if(gmain.p) return 1; + else return 0; +} + +void plumed_c2f(plumed p,char*c) { + unsigned i; + unsigned char* cc; + /* + Convert the address stored in p.p into a proper FORTRAN string + made of only ASCII characters. For this to work, the two following + assertions should be satisfied: + */ + assert(CHAR_BIT<=12); + assert(sizeof(p.p)<=16); + + assert(c); + cc=(unsigned char*)&p.p; + for(i=0; i. ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */ +#ifndef __PLUMED_wrapper_Plumed_h +#define __PLUMED_wrapper_Plumed_h + +/** +\page ReferencePlumedH Reference for interfacing MD codes with PLUMED + + Plumed.h and Plumed.c contain the external plumed interface, which is used to + integrate it with MD engines. This interface is very general, and is expected + not to change across plumed versions. Plumed.c also implements a dummy version + of the interface, so as to allow a code to be fully linked even if the plumed + library is not available yet. These files could be directly included in the official + host MD distribution. In this manner, it will be sufficient to link the plumed + library at link time (on all systems) or directly at runtime (on system where + dynamic loading is enabled) to include plumed features. + + Why is Plumed.c written in C and not C++? The reason is that the resulting Plumed.o + needs to be linked with the host MD code immediately (whereas the rest of plumed + could be linked a posteriori). Imagine the MD code is written in FORTRAN: when we + link the Plumed.o file we would like not to need any C++ library linked. In this + manner, we do not need to know which C++ compiler will be used to compile plumed. + The C++ library is only linked to the "rest" of plumed, which actually use it. + Anyway, Plumed.c is written in such a manner to allow its compilation also in C++ + (C++ is a bit stricter than C; compatibility is checked when PlumedStatic.cpp, + which basically includes Plumed.c, is compiled with the C++ compiler). This will + allow e.g. MD codes written in C++ to just incorporate Plumed.c (maybe renamed into + Plumed.cpp), without the need of configuring a plain C compiler. + + Plumed interface can be used from C, C++ and FORTRAN. Everything concerning plumed + is hidden inside a single object type, which is described in C by a structure + (struct \ref plumed), in C++ by a class (PLMD::Plumed) and in FORTRAN by a + fixed-length string (CHARACTER(LEN=32)). Obviously C++ can use both struct + and class interfaces, but the first should be preferred. The reference interface + is the C one, whereas FORTRAN and C++ interfaces are implemented as wrappers + around it. + + In the C++ interface, all the routines are implemented as methods of PLMD::Plumed. + In the C and FORTRAN interfaces, all the routines are named plumed_*, to + avoid potential name clashes. Notice that the entire plumed library + is implemented in C++, and it is hidden inside the PLMD namespace. + If the used C++ compiler supports C++11, PLMD::Plumed object defines move semantics + so as to be usable in STL containers. That is, you can declare a std::vector. + + Handlers to the plumed object can be converted among different representations, + to allow inter-operability among languages. In C, there are tools to convert + to/from FORTRAN, whereas in C++ there are tools to convert to/from FORTRAN and C. + + These handlers only contain a pointer to the real structure, so that + when a plumed object is brought from one language to another, + it brings a reference to the same environment. + + Moreover, to simplify life in all cases where a single Plumed object is + required for the entire simulation (which covers most of the practical + applications with conventional MD codes) it is possible to take advantage + of a global interface, which is implicitly referring to a unique global instance. + The global object should still be initialized and finalized properly. + + The basic method to send a message to plumed is +\verbatim + (C) plumed_cmd + (C++) PLMD::Plumed::cmd + (FORTRAN) PLUMED_F_CMD +\endverbatim + + To initialize a plumed object, use: +\verbatim + (C) plumed_create + (C++) (constructor of PLMD::Plumed) + (FORTRAN) PLUMED_F_CREATE +\endverbatim + + To finalize it, use +\verbatim + (C) plumed_finalize + (C++) (destructor of PLMD::Plumed) + (FORTRAN) PLUMED_F_FINALIZE +\endverbatim + + To access to the global-object, use +\verbatim + (C) plumed_gcreate, plumed_gfinalize, plumed_gcmd + (C++) PLMD::Plumed::gcreate, PLMD::Plumed::gfinalize, PLMD::Plumed::gcmd + (FORTRAN) PLUMED_F_GCREATE, PLUMED_F_GFINALIZE, PLUMED_F_GCMD +\endverbatim + + To check if the global object has been initialized, use +\verbatim + (C) plumed_ginitialized + (C++) PLMD::Plumed::ginitialized + (FORTRAN) PLUMED_F_GINITIALIZED +\endverbatim + + To check if plumed library is available (this is useful for runtime linking), use +\verbatim + (C) plumed_installed + (C++) PLMD::Plumed::installed + (FORTRAN) PLUMED_F_INSTALLED +\endverbatim + + To convert handlers use +\verbatim + (C) plumed_c2f (C to FORTRAN) + (C) plumed_f2c (FORTRAN to C) + (C++) Plumed(plumed) constructor (C to C++) + (C++) operator plumed() cast (C++ to C) + (C++) Plumed(char*) constructor (FORTRAN to C++) + (C++) toFortran(char*) (C++ to FORTRAN) +\endverbatim + +\verbatim + FORTRAN interface + SUBROUTINE PLUMED_F_INSTALLED(i) + INTEGER, INTENT(OUT) :: i + SUBROUTINE PLUMED_F_GINITIALIZED(i) + INTEGER, INTENT(OUT) :: i + SUBROUTINE PLUMED_F_GCREATE() + SUBROUTINE PLUMED_F_GCMD(key,val) + CHARACTER(LEN=*), INTENT(IN) :: key + UNSPECIFIED_TYPE, INTENT(INOUT) :: val(*) + SUBROUTINE PLUMED_F_GFINALIZE() + SUBROUTINE PLUMED_F_GLOBAL(p) + CHARACTER(LEN=32), INTENT(OUT) :: p + SUBROUTINE PLUMED_F_CREATE(p) + CHARACTER(LEN=32), INTENT(OUT) :: p + SUBROUTINE PLUMED_F_CMD(p,key,val) + CHARACTER(LEN=32), INTENT(IN) :: p + CHARACTER(LEN=*), INTENT(IN) :: key + UNSPECIFIED_TYPE, INTENT(INOUT) :: val(*) + SUBROUTINE PLUMED_F_FINALIZE(p) + CHARACTER(LEN=32), INTENT(IN) :: p +\endverbatim + + The main routine is "cmd", which accepts two arguments: + key is a string containing the name of the command + val is the argument. it is declared const so as to use allow passing const objects, but in practice plumed + is going to modify val in several cases (using a const_cast). + In some cases val can be omitted: just pass a NULL pointer (in C++, val is optional and can be omitted). + The set of possible keys is the real API of the plumed library, and will be expanded with time. + New commands will be added, but backward compatibility will be retained as long as possible. + + To pass plumed a callback function use the following syntax (not available in FORTRAN yet) +\verbatim + plumed_function_holder ff; + ff.p=your_function; + plumed_cmd(plumed,"xxxx",&ff); +\endverbatim + (this is passing the your_function() function to the "xxxx" command) +*/ + +#ifdef __cplusplus +extern "C" { +#endif + +/* Generic function pointer */ +typedef void (*plumed_function_pointer)(void); + +/** + \brief Holder for function pointer. + + To pass plumed a callback function use the following syntax: +\verbatim + plumed_function_holder ff; + ff.p=your_function; + plumed_cmd(plumed,"xxxx",&ff); +\endverbatim + (this is going to pass the your_function() function to the "xxxx" command) +*/ + +typedef struct { + plumed_function_pointer p; +} plumed_function_holder; + +/** + \brief Main plumed object + + This is an object containing a Plumed instance, which should be used in + the MD engine. It should first be initialized with plumed_create(), + then it communicates with the MD engine using plumed_cmd(). Finally, + before the termination, it should be deallocated with plumed_finalize(). + Its interface is very simple and general, and is expected + not to change across plumed versions. See \ref ReferencePlumedH. +*/ +typedef struct { + /** + \private + \brief Void pointer holding the real PlumedMain structure + */ + void*p; +} plumed; + +/** \relates plumed + \brief Constructor + + \return The constructed plumed object +*/ +plumed plumed_create(void); + +/** \relates plumed + \brief Tells p to execute a command + + \param p The plumed object on which command is acting + \param key The name of the command to be executed + \param val The argument. It is declared as const to allow calls like plumed_cmd(p,"A","B"), + but for some choice of key it can change the content +*/ +void plumed_cmd(plumed p,const char*key,const void*val); + +/** \relates plumed + \brief Destructor + + \param p The plumed object to be deallocated +*/ +void plumed_finalize(plumed p); + +/** \relates plumed + \brief Check if plumed is installed (for runtime binding) + + \return 1 if plumed is installed, 0 otherwise +*/ +int plumed_installed(void); + +/** \relates plumed + \brief Retrieves an handler to the global structure. +*/ +plumed plumed_global(void); + +/** \relates plumed + \brief Check if the global interface has been initialized + + \return 1 if plumed has been initialized, 0 otherwise +*/ +int plumed_ginitialized(void); + +/* global C interface, working on a global object */ + +/** \relates plumed + \brief Constructor for the global interface. + + \note Equivalent to plumed_create(), but initialize the static global plumed object +*/ +void plumed_gcreate(void); + +/** \relates plumed + \brief Tells to the global interface to execute a command. + + \param key The name of the command to be executed + \param val The argument. It is declared as const to allow calls like plumed_gcmd("A","B"), + but for some choice of key it can change the content + + \note Equivalent to plumed_cmd(), but acting on the global plumed object. + It thus does not require the plumed object to be specified. +*/ +void plumed_gcmd(const char* key,const void* val); + +/** \relates plumed + \brief Destructor for the global interface. + + \note Equivalent to plumed_finalize(), but acting on the global plumed object. + It thus does not require the plumed object to be specified. +*/ +void plumed_gfinalize(void); + +/* routines to convert char handler from/to plumed objects */ + +/** \related plumed + \brief Converts a C handler to a FORTRAN handler + + \param p The C handler + \param c The FORTRAN handler (a char[32]) + + This function can be used to convert a plumed object created in C to + a plumed handler that can be used in FORTRAN. +\verbatim +#include +int main(int argc,char*argv[]){ + plumed p; + p=plumed_create(); + char fortran_handler[32]; + plumed_c2f(p,fortran_handler); + printf("DEBUG: this is a string representation for the plumed handler: %s\n",fortran_handler); + fortran_routine(fortran_handler); + plumed_finalize(p); + return 0; +} +\endverbatim + Here `fortran_routine` is a routine implemented in FORTRAN that manipulates the + fortran_handler. +*/ +void plumed_c2f(plumed p,char* c); + +/** \related plumed + \brief Converts a FORTRAN handler to a C handler + \param c The FORTRAN handler (a char[32]) + \return The C handler + + This function can be used to convert a plumed object created in FORTRAN + to a plumed handler that can be used in C. +\verbatim +void c_routine(char handler[32]){ + plumed p; + p=plumed_f2c(handler); + plumed_cmd(p,"init",NULL); +} +\endverbatim + Here `c_routine` is a C function that can be called from FORTRAN + and interact with the provided plumed handler. +*/ +plumed plumed_f2c(const char* c); + +#ifdef __cplusplus +} +#endif + +#ifdef __cplusplus + +/* this is to include the NULL pointer */ +#include + +/* C++ interface is hidden in PLMD namespace (same as plumed library) */ +namespace PLMD { + +/** + C++ wrapper for \ref plumed. + + This class provides a C++ interface to PLUMED. +*/ + +class Plumed { + /** + C structure. + */ + plumed main; + /** + keeps track if the object was created from scratch using + the defaults destructor (reference=false) or if it was imported + from C or FORTRAN (reference=true). In the latter case, the + plumed_finalize() method is not called when destructing the object, + since it is expected to be finalized in the C/FORTRAN code + */ + bool reference; +public: + /** + Check if plumed is installed (for runtime binding) + \return true if plumed is installed, false otherwise + \note Equivalent to plumed_installed() but returns a bool + */ + static bool installed(); + /** + Check if global-plumed has been initialized + \return true if global plumed object (see global()) is initialized (i.e. if gcreate() has been + called), false otherwise. + \note Equivalent to plumed_ginitialized() but returns a bool + */ + static bool ginitialized(); + /** + Initialize global-plumed. + \note Equivalent to plumed_gcreate() + */ + static void gcreate(); + /** + Send a command to global-plumed + \param key The name of the command to be executed + \param val The argument. It is declared as const to allow calls like gcmd("A","B"), + but for some choice of key it can change the content + \note Equivalent to plumed_gcmd() + */ + static void gcmd(const char* key,const void* val); + /** + Finalize global-plumed + */ + static void gfinalize(); + /** + Returns the Plumed global object + \return The Plumed global object + */ + static Plumed global(); + /** + Constructor. + \note Performs the same task a plumed_create() + */ + Plumed(); + /** + Clone a Plumed object from a FORTRAN char* handler + \param c The FORTRAN handler (a char[32]). + + \attention The Plumed object created in this manner + will not finalize the corresponding plumed structure. + It is expected that the FORTRAN code calls plumed_c_finalize for it + */ +// to have maximum portability of this file I do not use the explicit keyword here +// I thus add a suppress command for cppcheck +// cppcheck-suppress noExplicitConstructor + Plumed(const char*c); + /** + Clone a Plumed object from a C plumed structure + \param p The C plumed structure. + + \attention The Plumed object created in this manner + will not finalize the corresponding plumed structure. + It is expected that the C code calls plumed_finalize for it + */ +// to have maximum portability of this file I do not use the explicit keyword here +// I thus add a suppress command for cppcheck +// cppcheck-suppress noExplicitConstructor + Plumed(plumed p); +private: + /** Copy constructor is disabled (private and unimplemented) + The problem here is that after copying it will not be clear who is + going to finalize the corresponding plumed structure. + */ + Plumed(const Plumed&); + /** Assignment operator is disabled (private and unimplemented) + The problem here is that after copying it will not be clear who is + going to finalize the corresponding plumed structure. + */ + Plumed&operator=(const Plumed&); +public: + /* + PLUMED 2.4 requires a C++11 compiler. + Anyway, since Plumed.h file might be redistributed with other codes + and it should be possible to combine it with earlier PLUMED versions, + we here explicitly check if C+11 is available before enabling move semantics. + This could still create problems if a compiler 'cheats', setting __cplusplus > 199711L + but not supporting move semantics. Hopefully will not happen! + */ +#if __cplusplus > 199711L + /** Move constructor. + Only if move semantics is enabled. + It allows storing PLMD::Plumed objects in STL containers. + */ + Plumed(Plumed&&); + /** Move assignment. + Only if move semantics is enabled. + */ + Plumed& operator=(Plumed&&); +#endif + /** + Retrieve the C plumed structure for this object + */ + operator plumed()const; + /** + Retrieve a FORTRAN handler for this object + \param c The FORTRAN handler (a char[32]). + */ + void toFortran(char*c)const; + /** + Send a command to this plumed object + \param key The name of the command to be executed + \param val The argument. It is declared as const to allow calls like p.cmd("A","B"), + but for some choice of key it can change the content + \note Equivalent to plumed_cmd() + */ + void cmd(const char*key,const void*val=NULL); + /** + Destructor + + Destructor is virtual so as to allow correct inheritance from Plumed object. + To avoid linking problems with g++, I specify "inline" also here (in principle + it should be enough to specify it down in the definition of the function, but + for some reason that I do not understand g++ does not inline it properly in that + case and complains when Plumed.h is included but Plumed.o is not linked. Anyway, the + way it is done here seems to work properly). + */ + inline virtual ~Plumed(); +}; + +/* All methods are inlined so as to avoid the compilation of an extra c++ file */ + +inline +bool Plumed::installed() { + return plumed_installed(); +} + +inline +Plumed::Plumed(): + main(plumed_create()), + reference(false) +{} + +inline +Plumed::Plumed(const char*c): + main(plumed_f2c(c)), + reference(true) +{} + +inline +Plumed::Plumed(plumed p): + main(p), + reference(true) +{} + +#if __cplusplus > 199711L +inline +Plumed::Plumed(Plumed&& p): + main(p.main), + reference(p.reference) +{} + +inline +Plumed& Plumed::operator=(Plumed&& p) { + main=p.main; + reference=p.reference; + return *this; +} +#endif + +inline +Plumed::operator plumed()const { + return main; +} + +inline +void Plumed::toFortran(char*c)const { + plumed_c2f(main,c); +} + +inline +void Plumed::cmd(const char*key,const void*val) { + plumed_cmd(main,key,val); +} + +inline +Plumed::~Plumed() { + if(!reference)plumed_finalize(main); +} + +inline +bool Plumed::ginitialized() { + return plumed_ginitialized(); +} + +inline +void Plumed::gcreate() { + plumed_gcreate(); +} + +inline +void Plumed::gcmd(const char* key,const void* val) { + plumed_gcmd(key,val); +} + +inline +void Plumed::gfinalize() { + plumed_gfinalize(); +} + +inline +Plumed Plumed::global() { + return plumed_global(); +} + +} + +#endif + + +#endif diff --git a/src/USER-PLUMED/Plumed.inc b/src/USER-PLUMED/Plumed.inc new file mode 100644 index 0000000000..26bc028d66 --- /dev/null +++ b/src/USER-PLUMED/Plumed.inc @@ -0,0 +1,3 @@ +# PLUMED: runtime installation +PLUMED_LOAD= /data/gt/mycodes/plumed2/src/lib/libplumedWrapper.a -ldl +PLUMED_DEPENDENCIES= diff --git a/src/USER-PLUMED/README b/src/USER-PLUMED/README new file mode 100644 index 0000000000..02a0229e57 --- /dev/null +++ b/src/USER-PLUMED/README @@ -0,0 +1,55 @@ +This package implements the "fix plumed" command, which can be used +in a LAMMPS input script. + +The fix allows enhanced sampling methods such as umbrella sampling and +metadynamics to be used. Furthermore, PLUMED can be used to perform a +wide range of analyses on trajectories on the fly as they are generated. + +The package uses the "PLUMED" library, whose source code is not included +in the LAMMPS source code distribution. The files in the USER-PLUMED package +folder implement an interface between LAMMPS and PLUMED, that are written +and maintained by Gareth Tribello (gareth.tribello@gmail.com). + +PLUMED must instead be downloaded and compiled separately to LAMMPS. This building +and compiling of PLUEMD can be done before or after the building of LAMMPS as LAMMPS +calls PLUMED as a dynamic library. If you wish to use PLUMED with LAMMPS, however, +you must run the command: + +make yes-user-plumed + +before compiling LAMMPS. Furthermore, you must ensure that PLUMED is in your +PATH when running a LAMMPS calculation that takes advantage of PLUMED. If +PLUMED is not in the PATH an error will be returned whenever LAMMPS encounters +the fix plumed command in its input. To be clear, however, LAMMPS will run if +it is compiled with fix-plumed enabled on inputs that do not contain a fix +plumed command when PLUMED is not in the PATH. + +More info about the PLUMED library can be found at: + +www.plumed.org + +and in the reference articles: + +PLUMED2: New feathers for an old bird +G.A. Tribello, M. Bonomi, D. Branduardi, C. Camilloni and G. Bussi, +Comp. Phys. Comm 185, 604 (2014) +https://doi.org/10.1016/j.cpc.2013.09.018 + +PLUMED: a portable plugin for free energy calculations with molecular dynamics +M. Bonomi, D. Branduardi, G. Bussi, C. Camilloni, D. Provasi, P. Raiteri, D. Donadio, F. Marinelli, F. Pietrucci, R.A. Broglia and M. Parrinello +Comp. Phys. Comm. 180, 1961 (2009) +https://doi.org/10.1016/j.cpc.2009.05.011 + +Instructions explaining how to use PLUMED and LAMMPS in tandem can be found on the PLUMED website, which also gives +numerous example scripts for PLUMED as well as citations to articles that dcoment the various methods that are +implemented within PLUMED. + +There are also example scripts for using this package in the folder +examples/USER/colvars, as well as on the GitHub page for PLUMED. + +Please contact Gareth Tribello (gareth.tribello@gmail.com) for questions +regarding this package. + +--------------------------------- + +Version: 2016-12-22 diff --git a/src/USER-PLUMED/fix_plumed.cpp b/src/USER-PLUMED/fix_plumed.cpp new file mode 100644 index 0000000000..91afd367b1 --- /dev/null +++ b/src/USER-PLUMED/fix_plumed.cpp @@ -0,0 +1,291 @@ +#include "math.h" +#include "stdlib.h" +#include "string.h" +#include "atom.h" +#include "update.h" +#include "force.h" +#include "respa.h" +#include "domain.h" +#include "error.h" +#include "group.h" +#include "fix_plumed.h" +#include "universe.h" +#include "compute.h" +#include "modify.h" +#include "pair.h" + +using namespace LAMMPS_NS; +using namespace PLMD; +using namespace FixConst; + +#define INVOKED_SCALAR 1 + +FixPlumed::FixPlumed(LAMMPS *lmp, int narg, char **arg) : + Fix(lmp, narg, arg), + p(NULL), + nlocal(0), + gatindex(NULL), + masses(NULL), + charges(NULL) +{ +// Not sure this is really necessary: + if (!atom->tag_enable) error->all(FLERR,"fix plumed requires atom tags"); +// Initialize plumed: + p=new PLMD::Plumed; + +// If the -partition option is activated then enable inter-partition communication + if (universe->existflag == 1) { + int me; + MPI_Comm inter_comm; + MPI_Comm_rank(world,&me); + // Change MPI_COMM_WORLD to universe->uworld which seems more appropriate + MPI_Comm_split(universe->uworld,me,0,&inter_comm); + p->cmd("GREX setMPIIntracomm",&world); + if (me == 0) { + // The inter-partition communicator is only defined for the root in + // each partition (a.k.a. world). This is due to the way in which + // it is defined inside plumed. + p->cmd("GREX setMPIIntercomm",&inter_comm); + } + p->cmd("GREX init",NULL); + } + // The general communicator is independent of the existence of partitions, + // if there are partitions, world is defined within each partition, + // whereas if partitions are not defined then world is equal to MPI_COMM_WORLD. + p->cmd("setMPIComm",&world); + +// Set up units +// LAMMPS units wrt kj/mol - nm - ps +// Set up units + + if (force->boltz == 1.0){ +// LAMMPS units lj + p->cmd("setNaturalUnits"); + } else { + double energyUnits=1.0; + double lengthUnits=1.0; + double timeUnits=1.0; + if (force->boltz == 0.0019872067){ +// LAMMPS units real :: kcal/mol; angstrom; fs + energyUnits=4.184; + lengthUnits=0.1; + timeUnits=0.001; + } else if (force->boltz == 8.617343e-5){ +// LAMMPS units metal :: eV; angstrom; ps + energyUnits=96.48530749925792; + lengthUnits=0.1; + timeUnits=1.0; + } else if (force->boltz == 1.3806504e-23){ +// LAMMPS units si :: Joule, m; s + energyUnits=0.001; + lengthUnits=1.e-9; + timeUnits=1.e-12; + } else if (force->boltz == 1.3806504e-16){ +// LAMMPS units cgs :: erg; cms;, s + energyUnits=6.0221418e13; + lengthUnits=1.e-7; + timeUnits=1.e-12; + } else if (force->boltz == 3.16681534e-6){ +// LAMMPS units electron :: Hartree, bohr, fs + energyUnits=2625.5257; + lengthUnits=0.052917725; + timeUnits=0.001; + } else error->all(FLERR,"Odd LAMMPS units, plumed cannot work with that"); + p->cmd("setMDEnergyUnits",&energyUnits); + p->cmd("setMDLengthUnits",&lengthUnits); + p->cmd("setMDTimeUnits",&timeUnits); + } + +// Read fix parameters: + int next=0; + for(int i=3;iexistflag == 1){ + // Each replica writes an independent log file + // with suffix equal to the replica id + char str_num[32], logFile[1024]; + sprintf(str_num,".%d",universe->iworld); + strncpy(logFile,arg[i],1024-32); + strcat(logFile,str_num); + p->cmd("setLogFile",logFile); + next=0; + } else { + // partition option not used + p->cmd("setLogFile",arg[i]); + next=0; + } + } + else if(!strcmp(arg[i],"plumedfile"))next=2; + else if(next==2){ + p->cmd("setPlumedDat",arg[i]); + next=0; + } + else error->all(FLERR,"syntax error in fix plumed - use 'fix name plumed plumedfile plumed.dat outfile plumed.out' "); + } + if(next==1) error->all(FLERR,"missing argument for outfile option"); + if(next==2) error->all(FLERR,"missing argument for plumedfile option"); + + p->cmd("setMDEngine","LAMMPS"); + + int natoms=int(atom->natoms); + p->cmd("setNatoms",&natoms); + + double dt=update->dt; + p->cmd("setTimestep",&dt); + + virial_flag=1; + scalar_flag = 1; + +// This is the real initialization: + p->cmd("init"); + +// Define compute to calculate potential energy + char *id_pe = (char *) "thermo_pe"; + int ipe = modify->find_compute(id_pe); + c_pe = modify->compute[ipe]; + // Trigger computation of potential energy every step + c_pe->addstep(update->ntimestep+1); +} + +FixPlumed::~FixPlumed() +{ + delete p; +} + +int FixPlumed::setmask() +{ + // set with a bitmask how and when apply the force from plumed + int mask = 0; + mask |= POST_FORCE; + mask |= THERMO_ENERGY; + mask |= POST_FORCE_RESPA; + mask |= MIN_POST_FORCE; + return mask; +} + +void FixPlumed::init() +{ + if (strcmp(update->integrate_style,"respa") == 0) + nlevels_respa = ((Respa *) update->integrate)->nlevels; +} + +void FixPlumed::setup(int vflag) +{ + if (strcmp(update->integrate_style,"verlet") == 0) + post_force(vflag); + else { + ((Respa *) update->integrate)->copy_flevel_f(nlevels_respa-1); + post_force_respa(vflag,nlevels_respa-1,0); + ((Respa *) update->integrate)->copy_f_flevel(nlevels_respa-1); + } +} + +void FixPlumed::min_setup(int vflag) +{ + post_force(vflag); +} + +void FixPlumed::post_force(int vflag) +{ + int update_gatindex=0; +// Try to find out if the domain decomposition has been updated: + if(nlocal!=atom->nlocal){ + if(charges) delete [] charges; + if(masses) delete [] masses; + if(gatindex) delete [] gatindex; + nlocal=atom->nlocal; + gatindex=new int [nlocal]; + masses=new double [nlocal]; + charges=new double [nlocal]; + update_gatindex=1; + } else { + for(int i=0;itag[i]-1){ + update_gatindex=1; + break; + } + } + } + MPI_Allreduce(MPI_IN_PLACE,&update_gatindex,1,MPI_INT,MPI_SUM,world); + +// In case it has been updated, rebuild the local mass/charges array +// and tell plumed about the change: + if(update_gatindex){ + for(int i=0;itag[i]-1; + masses[i]=atom->mass[atom->type[i]]; + if(atom->q) charges[i]=atom->q[atom->type[i]]; + } + p->cmd("setAtomsNlocal",&nlocal); + p->cmd("setAtomsGatindex",gatindex); + } + + +// set up local virial/box. plumed uses full 3x3 matrices + double virial[3][3]; + for(int i=0;i<3;i++) for(int j=0;j<3;j++) virial[i][j]=0.0; + double box[3][3]; + for(int i=0;i<3;i++) for(int j=0;j<3;j++) box[i][j]=0.0; + box[0][0]=domain->h[0]; + box[1][1]=domain->h[1]; + box[2][2]=domain->h[2]; + box[2][1]=domain->h[3]; + box[2][0]=domain->h[4]; + box[1][0]=domain->h[5]; + +// local variable with timestep: + int step=update->ntimestep; + +// pass all pointers to plumed: + p->cmd("setStep",&step); + p->cmd("setPositions",&atom->x[0][0]); + p->cmd("setBox",&box[0][0]); + p->cmd("setForces",&atom->f[0][0]); + p->cmd("setMasses",&masses[0]); + if(atom->q) p->cmd("setCharges",&charges[0]); + p->cmd("setVirial",&virial[0][0]); + p->cmd("getBias",&bias); + +// pass the energy + double pot_energy = 0.; + c_pe->compute_scalar(); + c_pe->invoked_flag |= INVOKED_SCALAR; + pot_energy = c_pe->scalar; + int nprocs; + // Divide energy by number of processors + // Plumed wants it this way + MPI_Comm_size(world,&nprocs); + pot_energy /= nprocs; + p->cmd("setEnergy",&pot_energy); + // Trigger computation of potential energy every step + c_pe->addstep(update->ntimestep+1); + +// do the real calculation: + p->cmd("calc"); + +// retransform virial to lammps representation: + Fix::virial[0]=-virial[0][0]; + Fix::virial[1]=-virial[1][1]; + Fix::virial[2]=-virial[2][2]; + Fix::virial[3]=-virial[0][1]; + Fix::virial[4]=-virial[0][2]; + Fix::virial[5]=-virial[1][2]; +} + +void FixPlumed::post_force_respa(int vflag, int ilevel, int iloop) +{ + if (ilevel == nlevels_respa-1) post_force(vflag); +} + +void FixPlumed::min_post_force(int vflag) +{ + post_force(vflag); +} + +double FixPlumed::compute_scalar() +{ + return bias; +} + + diff --git a/src/USER-PLUMED/fix_plumed.h b/src/USER-PLUMED/fix_plumed.h new file mode 100644 index 0000000000..055dc6f1d4 --- /dev/null +++ b/src/USER-PLUMED/fix_plumed.h @@ -0,0 +1,53 @@ + +#ifdef FIX_CLASS + +FixStyle(plumed,FixPlumed) + +#else + +#ifndef LMP_FIX_PLUMED_H +#define LMP_FIX_PLUMED_H + +#include "fix.h" +#include "compute.h" +// the plumed header that defines the class// +#include "../Plumed.h" + +namespace LAMMPS_NS { + +class FixPlumed : public Fix { + public: + FixPlumed(class LAMMPS *, int, char **); + ~FixPlumed(); + int setmask(); + void init(); + void setup(int); + void min_setup(int); + void post_force(int); + void post_force_respa(int, int, int); + void min_post_force(int); + double compute_scalar(); + + private: +// pointer to plumed object: + PLMD::Plumed*p; +// number of atoms local to this process: + int nlocal; +// array of atom indexes local to this process: + int*gatindex; +// array of masses for local atoms: + double*masses; +// array of charges for local atoms: + double*charges; +// this is something to enable respa + int nlevels_respa; +// output bias potential + double bias; +// Compute for the energy + class Compute *c_pe; +}; + +}; + +#endif +#endif