Restructure examples folder and add better READMEs (i-pi#245)

This is a first step towards better / more useful examples for i-PI.

We have restructured the folders and added lots of READMEs and documentation explaining what every example does. There is still lots to do so we'll keep the branch alive and re-open a new PR, but the very long list of comments below tells it's time to merge.

* Create new folder structure and move all examples to codes/ for the moment

* Quick-fix to the broken symlinks

* added the content of the PR i-pi#152

* restructure example folders

* move VASP to example/client folder

* move cav-md example

* update README demos

* fix broken sym link and rename some files

* add fhi-example

* delete pdb file

* tidy a bit

* update readme

* add Instanton and geop example

* update README

* move all vasp example files to the correct folder

* rename folder

* move demos folder, features/advanced_combinations

* install sdgml for testing

* Fix some example test

* fix more example test

* fix more example test

* change paths in exampletools.py

* lint

* move broken examples to broken folder

* move broken examples to broken folder

* update sym links zundel

* do not test examples in folder broken

* restore more sym links

* more fixes

* Disable in the workflow packages that are not being used

* move clients examples from temp to clients

* move clients examples from temp to clients

* reorganiz instanton examples

* reorganize instanton

* move reply to features

* Fixing a init file that was missing information (?)

* moved the NaCl example to BFGS. The previous example had incorrect units
parsed form the positions and cell and it also did not optimize the
cell.

* the nacl initial cell was absurd. Replaced with a more sane structure.

* Added LAMMPS paracetamol example

* moved and modified h2o molecules geop with lammps

* included the run script for paracetamol geop

* mentioned the new SD optimization example in the README

* Created directory for phonons feature. Included README and an example
with a run script.

* added normal mode analysis example for h2o-molecule

* renamed phonons to harmonic calculations as phonon quasiparticules are
not used for molecules

* updated the README of the harmonic calculations dir

* added example for an enmfd normalmode calculations for a water molecule

* included a new dir for IMF, VSCF, SCP

* populated directory for vibrational self consistent calculations

* added feature directory for KMC

* added feature directory for committee's

* typo

* included feature directory for quantum dynamics

* some edits to folder structure

* modified folder structure for ML committees

* modified TRPMD example folder name for consistency

* Fixed error in tex compilation due to a non-escaped %

* Updated fhi-aims tests

* Fixed and clarified problems

* Glitch

* Reverting a change and explaining issue. Need to address extra_type error in dummy driver.

* proofread readme. renamed loong folder name

* Cleaned up and reordered several README and example files

Fixed broken link in librascal example

* some more renaming

* Added folder with "demo" examples

* Small fixes to README

* Added a stub of a UDS run on HPC

* Bugfix for a weird output format change.

* Added a couple of READMEs

* Forgot to link water

* Making sure it works on another system

* added a multi-node example

* added open-paths feature

* added folder for PIMD

* added standard PIMD example

* added Cayley PIMD example

* added a piglet example

* added some more PIM examples

* added broken symlinks to committee and trpmd examples

* added README for classcial molecular dynamics

* added classical lammps example directories

* Diatomic Morse (for exercises)

* Actually had forgotten to commit last changes

* Just added command line options

* Added short-range qtip4pf potential for MTS demos

* Added documentation

* Stupid bug

* Trying to fix plumed+remd

* Fix array shapes to be compatible with new plumed typechecks

* Example of WTE metadynamics

* Added wte example to exclude list because we don't have plumed on the CI

... otherwise we wouldn't have had the bug in the 1st place....

* removed debug print statement

* Update documentation on tests and examples

The manual referred to nosetests and to a long-gone structure of the examples folder. I removed everything to avoid confusion (since the manual had become inconsistent with the README. 
We should also update  the manual when we finish i-pi#245

* Updated location of new WTE plumed example

* Removed leftover files

* Adding some dev instructions to README file

* improve lammps README

* extend FHI-aims readme

* remove repeated xyz files in cp2k and replace them by symlinks

* add symlink to librascal example

* correct main README.md in examples folder

* move broken to temp

* add Readme for instanton feature examples

* improve some readme files

* By general consensus, move demos outside examples

Also tried to update all references to the demos

* Some renaming, addressing issues in the review.

* Removed stray files and fixed links to zundel potential files

---------

Co-authored-by: Venkat Kapil <[email protected]>
Co-authored-by: Yair Litman <[email protected]>
Co-authored-by: Yair Litman <[email protected]>
Co-authored-by: Mariana Rossi <[email protected]>
Co-authored-by: Mariana Rossi 2 <[email protected]>

.github/workflows/examples.yml

@@ -39,6 +39,7 @@ jobs:
          pip install .
          pip install ase>=3.22.1
          echo "${GITHUB_WORKSPACE}/bin" >> $GITHUB_PATH
+         #pip install sgdml 
 
     - name: Install Gfortran and Compile Driver
       shell: bash

README.md

@@ -63,18 +63,23 @@ make
 cd ../..
 ```
 
-### Run one of the examples
+### Examples and demos
 
-This will first start the wrapper in the background, redirecting the output to
+The `examples` and `demos` folders contain inputs for many different types of
+calculations based on i-PI. Examples are typically minimal use-cases of specific
+features, while demos are more structured, tutorial-like examples that show how
+to realize more complex setups, and also provide a brief discussion of the 
+underlying algorithms.
+
+To run these examples, you should typically start i-PI, redirecting the output to
 a log file, and then run a couple of instances of the driver code. The progress
-of the wrapper is followed by monitoring the log file with the `tail` Linux
-command.
+of the wrapper is followed by monitoring the log file with the `tail` Linux command.
 
 Optionally, you can make a copy of the directory with the example somewhere
-else if you want to keep the i-PI directory clean.
+else if you want to keep the i-PI directory clean. For example
 
 ```bash
-cd examples/tutorial/tutorial-1/
+cd demos/para-h2-tutorial/tutorial-1/
 i-pi tutorial-1.xml > log &
 i-pi-driver -h localhost -p 31415 -m sg -o 15 &
 i-pi-driver -h localhost -p 31415 -m sg -o 15 &

demos/README.md

@@ -0,0 +1,30 @@
+i-PI demonstrations
+===================
+
+These examples provide more sophisticated and realistic simulations setups, that 
+demonstrate the usage of i-PI in a number of actual application scenarios. 
+This README provide an overview of the content of the various folders.
+
+
+p-H2: a tutorial
+----------------
+
+This is an introductory tutorial demonstrating classical and quantum simulations
+of para-hydrogen close to a phase transition. The example is discussed as a tutorial
+as part of the main documentations of i-PI. 
+
+
+Alchemical isotope exchanges
+----------------------------
+
+This example demonstrates the use of the alchemical isotope exchange method, which simplifies
+sampling equilibrium isotope fractionation between different molecular sites or different
+thermodynamic phases of matter.
+
+
+Kinetic Monte Carlo for Al6XXX
+------------------------------
+
+This is an experimental (and somewhat hacky) implementation of kinetic monte carlo, 
+specifically designed for the study of vacancy diffusion in aluminum alloys of the 6000 series (Al alloyed with ~1% Si and Mg).
+It demonstrates the use of i-PI with a ML potential, through LAMMPS patched with N2P2. Generalizing the implementation would require some coding, but most of the heavy lifting is already done. 

demos/al6xxx-kmc/README.md

@@ -0,0 +1,121 @@
+Kinetic Monte Carlo for vacancy diffusion in an Al-Si-Mg alloy
+==============================================================
+
+This example is somewhat atypical, in that it demonstrates a very ad hoc implementation of 
+KMC, designed specifically to run simulations of vacancy diffusion in an FCC alloy of aluminum
+with ~1% Si and Mg. 
+The example is interesting because it demonstrates several tricks that rely on the way i-PI treats
+atoms in communicating with the driver code, and uses a machine-learning potential to compute
+energy and forces. 
+The implmentation was first used in 
+[A. C. P. Jain, D. Marchand, A. Glensk, M. Ceriotti, and W. A. Curtin, "Machine learning for metallurgy III: A neural network potential for Al-Mg-Si," Phys. Rev. Materials 5(5), 053805 (2021)](https://doi.org/10.1103/PhysRevMaterials.5.053805), 
+and could probably be generalized with relatively little effort, so if you need to perform similar
+simulations for a different systems, read on. 
+
+Quasi-on-lattice KMC
+--------------------
+
+In general terms, KMC involves enumerating the possible reaction pathways for a given configuration
+of the system, estimating the associated rates, choosing one mechanism at random with the appropriate
+probability, and evolving the configuration accordingly. The process is repeated, giving a coarse-grained
+dynamics that can span much larger time scales than explicit molecular dynamics. 
+
+The implementation used here is peculiar, in that it considers an ideal fcc lattice decorated with
+atoms *and* vacancies, which define the configurations. However, the energies and barriers are estimated
+based on the energy of *relaxed* structures obtained starting from the ideal-lattice configuration.
+This approach allows incorporating the contributions from relaxation into the energetics of the KMC, while 
+still being able to enumerate possible reaction pathways based on the simple topology of the lattice. 
+
+Vacancies are assumed to be able to diffuse in all directions, and the barrier is taken to be a fixed 
+shift on top of the mean between reactant and product energies. 
+Thus, at each step the code has to estimate the energy of all the possible outcomes for vacancy 
+diffusion, and compute their relaxed energies. 
+
+
+Required software
+-----------------
+
+This simulation requires the use of [LAMMPS](https://www.lammps.org/), including the 
+`fix_ipi`, as well as compiling the [N2P2](https://compphysvienna.github.io/n2p2) library
+and its LAMMPS interface, to compute the Behler-Parrinello potential for Al-Mg-Si alloys. 
+
+Usage
+-----
+
+After compiling LAMMPS with the N2P2 module, You will have to launch i-PI
+
+```bash
+i-pi input.xml &> log.ipi
+```
+
+and then run one or more instances of LAMMPS
+
+```bash
+$LAMMPS_EXECUTABLE < in.lmp &> log.lammps
+```
+
+What to look for
+----------------
+
+The `N2P2` parameters (and the network weights) are in the folder `potential_n2p2_v2ag`. 
+The file `input.nn` contains all the symmetry-function parameters, as well as the parameters that 
+were used for training.
+
+
+The `motion` class that implements KMC is very verbose, and reflects the highly ad-hoc nature 
+of the implementation. It is important to understand that i-PI in this case will essentially ignore
+the initialization structure, and create a new one, with a `ncell×ncell×ncell` replication of a primitive 
+fcc cell, with the stated (conventional cell) lattice parameter `a0`. The cell defaults to a pure Al
+matrix, and the first atoms are set to be `nsi` Si atoms and `nmg` Mg atoms. The last `nvac` atoms are 
+set to vacancies (V), and we will discuss later how these are treated. 
+
+Then, the motion class specifies what geometry optimization algorithm should be run (these can be any 
+of the `optimizer` classes uses with `motion mode='geop'`) and how many steps of relaxation (`nstep`) should be performed after each KMC reaction. It then defines the estimated diffusion barrier and prefactor.
+
+Finally, the class implements some optimizations. First, `<neval>` potential reactive pathways are
+evaluated in parallel (provided enough calculators are available). Second, reactions are cached, 
+i.e. the relaxed positions and energies are saved and re-used in case the system visits again the same
+configuration. This is useful when the system gets trapped into a loop of states connected by low 
+reaction barriers. 
+
+```xml
+<motion mode='al-kmc'>
+    <al6xxx_kmc>
+       <ncell> 5 </ncell>
+       <a0 units='angstrom'> 4.057 </a0>       
+       <nsi> 6 </nsi>
+       <nmg> 6 </nmg>
+       <nvac> 2 </nvac>
+       <optimizer mode='lbfgs'> 
+            <ls_options> <iter> 3 </iter> </ls_options>
+       </optimizer>
+       <nstep> 5 </nstep>
+       <diffusion_barrier_al units='electronvolt'> 0.52  </diffusion_barrier_al> <!-- Mantina2009 dHm -->
+       <diffusion_prefactor_al units='terahertz'> 16.6 </diffusion_prefactor_al> <!-- Mantina2009 v*  -->
+       <neval> 8 </neval>
+       <ecache_file> KMC_ECACHE </ecache_file>
+       <qcache_file> KMC_QCACHE </qcache_file>
+       <max_cache_len> 100 </max_cache_len>   
+   </al6xxx_kmc>
+</motion>
+```
+
+This application is also a good example of how the communication between i-PI and the driver code can
+be exploited to realize complicated simulations. Here, we want to define some "virtual atoms" that 
+act as placeholders for vacancies, but LAMMPS (at least with the N2P2 plugin) has no mechanism to 
+interpret these as anything but actual atoms. We get around this problem by using the 
+`activelist` option to the forcefield command, that specifies the list of atoms that are
+actually passed on to the driver for calculation. So, even if the fcc box has 125 atoms,
+and on the i-PI side the last to 'V' atoms get printed out and manipulated, only the first
+123 "real" atoms are passed to LAMMPS. 
+
+```xml
+<activelist> [
+       0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,
+33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,56,57,58,59,60,61,62,63,64,
+65,66,67,68,69,70,71,72,73,74,75,76,77,78,79,80,81,82,83,84,85,86,87,88,89,90,91,92,93,94,95,96,
+97,98,99,100,101,102,103,104,105,106,107,108,109,110,111,112,113,114,115,116,117,118,119,120,121,122] </activelist>
+```
+
+Incidentally, this means that on the LAMMPS side we should have only 123 atoms: 
+6 Si, 6 Mg and then 111 Al atoms. 

examples/lammps/kmc-al6xxx/input.xml → demos/al6xxx-kmc/input.xml

@@ -1,47 +1,47 @@
-<!--- To run with a reasonable LJ potential use
+<!--- To run with a reasonable LJ potential for testing purposes, use
   i-pi-driver -m lj -u -h al6xxx -o 5.0,0.1,12.0
 -->
-<simulation verbosity="low">
-    <output prefix="simulation">
-        <properties stride="1" filename="out"> [ step, time{picosecond}, potential] </properties>
-        <properties stride="1" filename="vac"> [ step, time{picosecond}, atom_x{angstrom}(123), atom_x{angstrom}(124) ] </properties>
-        <trajectory filename="pos" stride="1" cell_units="angstrom" format="xyz" bead="0"> positions{angstrom} </trajectory>
+<simulation verbosity='low'>
+    <output prefix='simulation'>
+        <properties stride='1' filename='out'> [ step, time{picosecond}, potential] </properties>
+        <properties stride='1' filename='vac'> [ step, time{picosecond}, atom_x{angstrom}(123), atom_x{angstrom}(124) ] </properties>
+        <trajectory filename='pos' stride='1' cell_units='angstrom' format='xyz' bead='0'> positions{angstrom} </trajectory>
     </output>
     <total_steps> 2 </total_steps>
     <prng> <seed> 12345 </seed> </prng>
-    <ffsocket name="lmpserial" mode="unix">
+    <ffsocket name='lmpserial' mode='unix'>
         <address> localhost_0_0 </address> <latency> 1e-3 </latency>
         <activelist> [0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,56,57,58,59,60,61,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76,77,78,79,80,81,82,83,84,85,86,87,88,89,90,91,92,93,94,95,96,97,98,99,100,101,102,103,104,105,106,107,108,109,110,111,112,113,114,115,116,117,118,119,120,121,122] </activelist>
     </ffsocket>
     <system>
-        <initialize nbeads="1">
-            <file mode="xyz" units="angstrom"> init.xyz </file>
-            <velocities mode="thermal" units="kelvin"> 0 </velocities>
+        <initialize nbeads='1'>
+            <file mode='xyz' units='angstrom'> init.xyz </file>
+            <velocities mode='thermal' units='kelvin'> 0 </velocities>
         </initialize>
         <forces>
-            <force forcefield="lmpserial"> </force>
+            <force forcefield='lmpserial'> </force>
         </forces>
-        <motion mode="al-kmc">
+        <motion mode='al-kmc'>
             <al6xxx_kmc>
-               <optimizer mode="lbfgs"> 
+               <optimizer mode='lbfgs'> 
                     <ls_options> <iter> 3 </iter> </ls_options>
                </optimizer>
-               <a0 units="angstrom"> 4.057 </a0>
+               <a0 units='angstrom'> 4.057 </a0>
                <nstep> 5 </nstep>
                <ncell> 5 </ncell>
                <nsi> 6 </nsi>
                <nmg> 6 </nmg>
                <nvac> 2 </nvac>
                <neval> 8 </neval>
-	           <diffusion_barrier_al units="electronvolt"> 0.52  </diffusion_barrier_al> <!-- Mantina2009 dHm -->
-               <diffusion_prefactor_al units="terahertz"> 16.6 </diffusion_prefactor_al> <!-- Mantina2009 v*  -->
+	           <diffusion_barrier_al units='electronvolt'> 0.52  </diffusion_barrier_al> <!-- Mantina2009 dHm -->
+               <diffusion_prefactor_al units='terahertz'> 16.6 </diffusion_prefactor_al> <!-- Mantina2009 v*  -->
 	           <ecache_file> KMC_ECACHE </ecache_file>
                <qcache_file> KMC_QCACHE </qcache_file>
                <max_cache_len> 100 </max_cache_len>   
            </al6xxx_kmc>    
         </motion>
         <ensemble>
-            <temperature units="kelvin">300</temperature>
+            <temperature units='kelvin'>300</temperature>
         </ensemble>
     </system>
 </simulation>

examples/lammps/kmc-al6xxx/potential_n2p2_v2ag/input.nn → demos/al6xxx-kmc/potential_n2p2_v2ag/input.nn

@@ -6,8 +6,6 @@ elements                        Al Mg Si           # Specification of elements.
 atom_energy  Al    -19.6286626  #-2.4092354 
 atom_energy  Mg    -16.7493346  #-62.6068620
 atom_energy  Si    -5.5274864   #-5.5597835
-#atom_energy                     S  0.0         # Free atom reference energy (S).
-#atom_energy                     Cu 0.0         # Free atom reference energy (Cu).
 cutoff_type                     2          # Cutoff type (optional argument: shift parameter alpha).
 scale_symmetry_functions                       # Scale all symmetry functions with min/max values.
 #scale_symmetry_functions_sigma                 # Scale all symmetry functions with sigma.

demos/al6xxx-kmc/run.sh

@@ -0,0 +1,16 @@
+ipi=i-pi
+lmp=lmp_mpi 
+sleep_time=10
+
+${ipi} input.xml > log.i-pi & 
+echo "# i-PI is running"
+
+echo "# Waiting for ${sleep_time} (s) before executing LAMMPS"
+sleep ${sleep_time}
+
+${lmp} < in.lmp > /dev/null & 
+echo "# LAMMPS is running"
+
+wait
+
+echo "# Simulation complete"

demos/alchemical-isotope-exchanges/README.md

@@ -0,0 +1,72 @@
+Alchemical isotope exchanges
+============================
+
+This example demonstrates the use of the alchemical exchange method implemented in i-PI
+(cf. the theory introduced in [Jian Liu, Richard S Andino, Christina M Miller, Xin Chen, David M Wilkins, Michele Ceriotti, David E Manolopoulos, “A surface-specific isotope effect in mixtures of light and heavy water”, J. Phys. Chem. C 117(6), 2944-2951 (2013)](http://dx.doi.org/10.1021/jp311986m) [bibtex](https://www.doi2bib.org/bib/10.1021/jp311986m), 
+and the implementation discussed in [Cheng, Bingqing, Jörg Behler, Michele Ceriotti, “Nuclear Quantum Effects in Water at the Triple Point: Using Theory as a Link Between Experiments.” J. Phys. Chem. Lett. 7(12), 2210-2215 (2016)](https://ipi-code.org/about/features/dx.doi.org/10.1021/acs.jpclett.6b00729) [bibtex](http://www.doi2bib.org/bib/10.1021%2Facs.jpclett.6b00729) ).
+
+The masses of different isotopes (and the corresponding labels) are exchanged with a 
+Monte Carlo acceptance criterion based on the full path integral Hamiltonian, leading to
+equilibration between different molecular positions and/or thermodynamic phases. 
+
+In this example, the simulation starts with 32 water HOD molecules, and rapidly reaches
+an equilibrium with 1:2:1 H2O:HOD:D2O ratio; simulations at a lower temperature (and 
+using a higher number of PI replicas) will lead to a small deviation from this ratio,
+due to the presence of equilibrium isotope effects. 
+
+
+Required software
+-----------------
+
+This simulation requires the use of [LAMMPS](https://www.lammps.org/), including the 
+`fix_ipi` module and all the modules needed to use a TIP4P water potential.
+
+
+Usage
+-----
+
+Much as with any PIMD simulation, you will have to launch i-PI
+
+```bash
+i-pi input.xml &> log.ipi
+```
+
+and then run one or more instances of the driver code
+
+```bash
+$LAMMPS_EXECUTABLE < in.lmp &> log.lammps
+```
+
+Note that the simulation uses some common files from the `init_files` folder,
+so you may have to copy those and update paths if you want to run it outside
+the examples folder.
+
+
+What to look for
+----------------
+
+The `input.xml` file describes a run-of-the-mill PIMD calculation. The key section for
+this example is 
+
+```xml
+<motion mode='alchemy'>
+  <alchemy>
+    <names> [ H, D ] </names>
+    <nxc> 10 </nxc>
+  </alchemy>
+</motion>
+```
+
+That requires attempting on average 10 exhanges between H and D atoms per step. Given
+that the exchange is assumed to only change the mass matrix, the acceptance probability
+can be computed quickly, without having to re-evaluate the potential energy. 
+
+On the other hand, note that this assumes that the atoms being swapped are bona fide 
+isotopes, and that they behave chemically in the same way: in this example, the driver
+code does not need to know about the swap having happened, because H and D behave the same 
+from the point of view of the forcefield.
+
+If you are looking for actual chemical swaps, these can be realized using a different
+class, [atomswap](https://ipi-code.org/i-pi/input-reference.html#atomswap).
+
+See in the output positions file how the atomic labels are swapped between structures.

demos/alchemical-isotope-exchanges/in.lmp

@@ -0,0 +1,19 @@
+units		electron
+atom_style	full
+
+# q-TIP4P-f parameters
+pair_style      lj/cut/tip4p/long 1 2 1 1 0.278072379 17.007
+bond_style      class2 
+angle_style     harmonic
+kspace_style	pppm/tip4p 0.0001
+
+read_data	water_32_data.lmp
+pair_coeff  * * 0 0
+pair_coeff  1  1  0.000295147 5.96946
+
+neighbor	2.0 bin
+
+timestep	0.00025 # timestep is irrelevant: MD happens within i-PI
+fix 1 all ipi h2o-alchemical 32345 unix
+run		100000000
+