smooth version of deepmd. better mixed systems. add mpi support in lammps. add small box support.

Author: Han Wang

README.md

@@ -17,7 +17,9 @@
 - [License](#license)
 
 # Install DeePMD-kit
-The installation of the DeePMD-kit is lengthy, but do not be panic. Just follow step by step. Wish you good luck..
+The installation of the DeePMD-kit is lengthy, but do not be panic. Just follow step by step. Wish you good luck.. 
+
+A docker for installing the DeePMD-kit on CentOS 7 is available [here](https://github.com/TimChen314/deepmd-kit_docker).
 
 ## Install tensorflow's Python interface 
 There are two ways of installing the Python interface of tensorflow, either [using google's binary](https://www.tensorflow.org/install/install_linux), or [installing from sources](https://www.tensorflow.org/install/install_sources). When you are using google's binary, do not forget to add the option `-DTF_GOOGLE_BIN=true` when building DeePMD-kit.
@@ -168,12 +170,12 @@ dp_frz  dp_ipi  dp_mdnn  dp_test  dp_train
 ```
 
 ## Install Lammps' DeePMD-kit module
-DeePMD-kit provide module for running serial MD simulation with Lammps. Notice that the parallel running is not support at this moment. Now make the DeePMD-kit module for lammps.
+DeePMD-kit provide module for running MD simulation with Lammps. Now make the DeePMD-kit module for lammps.
 ```bash
 cd $deepmd_source_dir/source/build
 make lammps
 ```
-If everything works fine, DeePMD-kit will generate a module called `USER-DEEPMD` in the `build` directory. Now download your favorite Lammps code, and uncompress it (I assume that you have downloaded the tar `lammps-stable.tar.gz`)
+DeePMD-kit will generate a module called `USER-DEEPMD` in the `build` directory. Now download your favorite Lammps code, and uncompress it (I assume that you have downloaded the tar `lammps-stable.tar.gz`)
 ```bash
 cd /some/workspace
 tar xf lammps-stable.tar.gz
@@ -186,10 +188,16 @@ cp -r $deepmd_source_dir/source/build/USER-DEEPMD .
 Now build Lammps
 ```bash
 make yes-user-deepmd
-make serial -j4
+make mpi -j4
+```
+The option `-j4` means using 4 processes in parallel. You may want to be use a different number according to your hardware. 
+
+If everything works fine, you will end up with an executable `lmp_mpi`.
+
+The DeePMD-kit module can be removed from Lammps source code by 
+```bash
+make no-user-deepmd
 ```
-The option `-j4` means using 4 processes in parallel. You may want to be use a different number according to your hardware. If everything works fine, you will end up with an executable
-`lmp_serial`.
 
 # Use DeePMD-kit
 In this text, we will call the deep neural network that is used to represent the interatomic interactions (Deep Potential) the **model**. The typical procedure of using DeePMD-kit is 
@@ -238,6 +246,7 @@ box.raw  coord.raw  energy.raw  force.raw  set.000  set.001  set.002  type.raw
 It generates two sets `set.000`, `set.001` and `set.002`, with each set contains 2000 frames. The last set (`set.002`) is used as testing set, while the rest sets (`set.000` and `set.001`) are used as training sets. One do not need to take care the binary data files in each of the `set.*` directories. The path containing `set.*` and `type.raw` is called a *system*. 
 
 ## Train a model
+### The standard DeePMD model
 The method of training is explained in our [DeePMD paper][1]. With the source code we provide a small training dataset taken from 400 frames generated by NVT ab-initio water MD trajectory with 300 frames for training and 100 for testing. [An example training parameter file](./examples/train/water.json) is provided. One can try with the training by
 ```bash
 $ cd $deepmd_source_dir/examples/train/
@@ -247,10 +256,10 @@ $ $deepmd_root/bin/dp_train water.json
 ```json
 {
     "_comment": " model parameters",
+    "use_smooth":	false,
     "sel_a":		[16, 32],
     "sel_r":		[30, 60],
-    "rcut_a":		-1,
-    "rcut_r":		6.00,
+    "rcut":		6.00,
     "axis_rule":	[0, 1, 0, 0, 1, 1, 0, 0, 0, 0, 1, 0],
     "_comment":	" default rule: []",
     "_comment":	" user defined rule: for each type provides two axes, ",
@@ -273,10 +282,9 @@ $ $deepmd_root/bin/dp_train water.json
     "limit_pref_e":	8,
     "start_pref_f":	1000,
     "limit_pref_f":	1,
-    "start_pref_v":	0.02,
-    "limit_pref_v":	8,
+    "start_pref_v":	0,
+    "limit_pref_v":	0,
 
-    "num_threads":	4,
     "seed":		1,
 
     "_comment": " display and restart",
@@ -286,7 +294,6 @@ $ $deepmd_root/bin/dp_train water.json
     "numb_test":	100,
     "save_freq":	100,
     "save_ckpt":	"model.ckpt",
-    "restart":		false,
     "load_ckpt":	"model.ckpt",
     "disp_training":	true,
     "time_training":	true,
@@ -295,7 +302,7 @@ $ $deepmd_root/bin/dp_train water.json
 }
 ```
 
-The option **`rcut_r`** is the cut-off radius for neighbor searching. The `sel_a` and `sel_r` are the maximum selected numbers of fully-local-coordinate and radial-only-coordinate atoms from the neighbor list, respectively. `sel_a + sel_r` should larger than the maximum possible number of neighbors in the cut-off radius. `sel_a` and `sel_r` are vectors, the length of the vectors are same as the number of atom types in the system. `sel_a[i]` and `sel_r[i]` denote the selected number of neighbors of type `i`.
+The option **`rcut`** is the cut-off radius for neighbor searching. The `sel_a` and `sel_r` are the maximum selected numbers of fully-local-coordinate and radial-only-coordinate atoms from the neighbor list, respectively. `sel_a + sel_r` should larger than the maximum possible number of neighbors in the cut-off radius. `sel_a` and `sel_r` are vectors, the length of the vectors are same as the number of atom types in the system. `sel_a[i]` and `sel_r[i]` denote the selected number of neighbors of type `i`.
 
 The option **`axis_rule`** specifies how to make the axis for the local coordinate of each atom. For each atom type, 6 integers should be provided. The first three for the first axis, while the last three for the second axis. Within the three integers, the first one specifies if the axis atom is fully-local-coordinated (`0`) or radial-only-coordinated (`1`). The second integer specifies the type of the axis atom. If this number is less than 0, saying `t < 0`, then this axis exclude atom of type `-(t+1)`. If the third integer is, saying `s`, then the axis atom is the `s`th nearest neighbor satisfying the previous two conditions. 
 
@@ -314,15 +321,63 @@ The options **`start_pref_e`**, **`limit_pref_e`**, **`start_pref_f`**, **`limit
 ```math
 w_f(t) = start_pref_f * ( lr(t) / start_lr ) + limit_pref_f * ( 1 - lr(t) / start_lr )
 ```
+Since we do not have virial data, the virial prefactors `start_pref_v` and `limit_pref_v` are set to 0.
 
-The option **`num_threads`** specifies the number of threads used in the training. 
-
-The option **`seed`** specifies the random seed for neural network initialization. 
+The option **`seed`** specifies the random seed for neural network initialization. If not provided, the `seed` will be initialized with `None`.
 
 During the training, the error of the model is tested every **`disp_freq`** batches with **`numb_test`** frames from the last set in the **`systems`** directory on the fly, and the results are output to **`disp_file`**. 
 
 Checkpoints will be written to files with prefix **`save_ckpt`** every **`save_freq`** batches. If **`restart`** is set to `true`, then the training will start from the checkpoint named **`load_ckpt`**, rather than from scratch.
 
+Several command line options can be passed to `dp_train`, this can be checked with
+```bash
+$ $deepmd_root/bin/dp_train --help
+```
+An explanation will be provided
+```
+positional arguments:
+  INPUT                 the input json database
+
+optional arguments:
+  -h, --help            show this help message and exit
+  -t INTER_THREADS, --inter-threads INTER_THREADS
+                        With default value 0. Setting the "inter_op_parallelism_threads" key for the tensorflow, the "intra_op_parallelism_threads" will be set by the env variable OMP_NUM_THREADS
+  --init-model INIT_MODEL
+                        Initialize a model by the provided checkpoint
+  --restart RESTART     Restart the training from the provided checkpoint
+```
+The keys `intra_op_parallelism_threads` and `inter_op_parallelism_threads` are Tensorflow configurations for multithreading, which are explained [here](https://www.tensorflow.org/performance/performance_guide#optimizing_for_cpu). Skipping `-t` and `OMP_NUM_THREADS` leads to the default setting of these keys in the Tensorflow.
+
+**`--init-model model.ckpt`**, for example, initializes the model training with an existing model that is stored in the checkpoint `model.ckpt`, the network architectures should match.
+
+**`--restart model.ckpt`**, continues the training from the checkpoint `model.ckpt`.
+
+### The smooth DeePMD model
+The smooth version of DeePMD can be trained by the DeePMD-kit. [An example training parameter file](./examples/train/water_smth.json) is provided. One can try with the training by
+```bash
+$ cd $deepmd_source_dir/examples/train/
+$ $deepmd_root/bin/dp_train water_smth.json
+```
+The difference between the standard and smooth DeePMD models lies in the model parameters:
+```json
+    "use_smooth":	true,
+    "sel_a":		[46, 92],
+    "rcut_smth":	5.80,
+    "rcut":		6.00,
+    "filter_neuron":	[25, 50, 100],
+    "filter_resnet_dt":	false,
+    "n_axis_neuron":	16,
+    "n_neuron":		[240, 240, 240],
+    "resnet_dt":	true,
+```
+The `sel_r` option is skipped by the smooth version and the model use fully-local-coordinate for all neighboring atoms. The `sel_a` should larger than the maximum possible number of neighbors in the cut-off radius `rcut`. 
+
+The descriptors will decay smoothly from **`rcut_smth`** to the cutoff radius `rcut`.
+
+**`filter_neuron`** provides the size of the filter network (also called local-embedding network). If the size of the next layer is the same or twice as the previous layer, then a skip connection is build (ResNet). **`filter_resnet_dt`** tells if a timestep is used in the skip connection. By default it is `false`. **`n_axis_neuron`** specifies the number of axis filter, which should be much smaller than the size of the last layer of the filter network.
+
+**`n_neuron`** specifies the fitting network. If the size of the next layer is the same as the previous layer, then a skip connection is build (ResNet). **`resnet_dt`** tells if a timestep is used in the skip connection. By default it is `true`. 
+
 
 ## Freeze the model
 The trained neural network is extracted from a checkpoint and dumped into a database. This process is called "freeze" a model. Typically one does
@@ -331,7 +386,6 @@ $ $deepmd_root/bin/dp_frz -o graph.pb
 ```
 in the folder where the model is trained. The output database is called `graph.pb`.
 
-
 ## Run MD with Lammps
 Run an MD simulation with Lammps is simpler. In the Lammps input file, one needs to specify the pair style as follows
 ```bash
@@ -340,6 +394,17 @@ pair_coeff
 ```
 where `graph.pb` is the file name of the frozen model. The `pair_coeff` should be left blank. It should be noted that Lammps counts atom types starting from 1, therefore, all Lammps atom type will be firstly subtracted by 1, and then passed into the DeePMD-kit engine to compute the interactions.
 
+### With long-range interaction
+The reciprocal space part of the long-range interaction can be calculated by lammps command `kspace_style`. To use it with DeePMD-kit, one writes 
+```bash
+pair_style	hybrid/overlay deepmd graph.pb coul/long 9.0
+pair_coeff	* * deepmd
+pair_coeff	* * coul/long
+pair_modify	pair coul/long compute no
+kspace_style	pppm 1.0e-5
+kspace_modify	gewald 0.45
+```
+In this setting, the direct space part of the long-range interaction is ignored by the `pair_modify` command, because this part is fitted in the DeePMD model. The splitting parameter `gewald` is modified by the `kspace_modify` command.
 
 ## Run path-integral MD with i-PI
 The i-PI works in a client-server model. The i-PI provides the server for integrating the replica positions of atoms, while the DeePMD-kit provides a client named `dp_ipi` that computes the interactions (including energy, force and virial). The server and client communicates via the Unix domain socket or the Internet socket. The client can be started by

data/raw/copy_raw.py

@@ -7,9 +7,9 @@
 def copy (in_dir,
           out_dir,
           ncopies = [1,1,1]) :
-    has_energy = os.path.isfile (in_dir + "energy.raw")
-    has_force  = os.path.isfile (in_dir + "force.raw")
-    has_virial = os.path.isfile (in_dir + "virial.raw")
+    has_energy = os.path.isfile (in_dir + "/energy.raw")
+    has_force  = os.path.isfile (in_dir + "/force.raw")
+    has_virial = os.path.isfile (in_dir + "/virial.raw")
 
     i_box       = np.loadtxt (in_dir + "/box.raw")
     i_coord     = np.loadtxt (in_dir + "/coord.raw")
@@ -65,7 +65,7 @@ def copy (in_dir,
         np.savetxt (out_dir + "/force.raw",         o_force)
     if has_virial :
         np.savetxt (out_dir + "/virial.raw",        o_virial)
-    np.savetxt (out_dir + "/type.raw",          o_type)
+    np.savetxt (out_dir + "/type.raw",          o_type, fmt = '%d')
     np.savetxt (out_dir + "/ncopies.raw",       ncopies, fmt = "%d")
 
 def _main () :

data/raw/raw_to_set.sh

@@ -6,6 +6,7 @@ if test $# -ge 1; then
     nline_per_set=$1
 fi
 
+rm -fr set.*
 echo nframe is `cat energy.raw | wc -l`
 echo nline per set is $nline_per_set
 

data/raw/shuffle_raw.py

@@ -42,17 +42,25 @@ def _main () :
         print ("# no file to shuffle, exit")
         return
 
+    assert ("box.raw" in raws)
+    tmp = np.loadtxt(os.path.join(inpath, "box.raw"))
+    tmp = np.reshape(tmp, [-1, 9])
+    nframe = tmp.shape[0]
+    print(nframe)
+
     print ("# will shuffle raw files " + str(raws) + 
            " in dir " + inpath +
            " and output to dir " + outpath)
 
     tmp = np.loadtxt (inpath + "/" + raws[0])
+    tmp = np.reshape(tmp, [nframe, -1])
     nframe = tmp.shape[0]
     idx = np.arange (nframe)
     np.random.shuffle(idx)
 
     for ii in raws : 
         data = np.loadtxt(inpath + "/" + ii)
+        data = np.reshape(data, [nframe, -1])
         data = data [idx]
         np.savetxt (outpath + "/" + ii, data)
 

examples/train/water.json

@@ -1,9 +1,9 @@
 {
     "_comment": " model parameters",
+    "use_smooth":	false,
     "sel_a":		[16, 32],
     "sel_r":		[30, 60],
-    "rcut_a":		-1,
-    "rcut_r":		6.00,
+    "rcut":		6.00,
     "axis_rule":	[0, 1, 0, 0, 1, 1, 0, 0, 0, 0, 1, 0],
     "_comment":	" default rule: []",
     "_comment":	" user defined rule: for each type provides two axes, ",
@@ -29,7 +29,6 @@
     "start_pref_v":	0,
     "limit_pref_v":	0,
 
-    "num_threads":	4,
     "seed":		1,
 
     "_comment": " display and restart",
@@ -39,7 +38,6 @@
     "numb_test":	100,
     "save_freq":	100,
     "save_ckpt":	"model.ckpt",
-    "restart":		false,
     "load_ckpt":	"model.ckpt",
     "disp_training":	true,
     "time_training":	true,

examples/train/water_smth.json

@@ -0,0 +1,46 @@
+{
+    "_comment": " model parameters",
+    "use_smooth":	true,
+    "sel_a":		[46, 92],
+    "rcut_smth":	5.80,
+    "rcut":		6.00,
+    "filter_neuron":	[25, 50, 100],
+    "filter_resnet_dt":	false,
+    "n_axis_neuron":	16,
+    "n_neuron":		[240, 240, 240],
+    "resnet_dt":	true,
+
+    "_comment": " traing controls",
+    "systems":		["../data/water/"],
+    "set_prefix":	"set",    
+    "stop_batch":	1000000,
+    "batch_size":	1,
+    "start_lr":		0.005,
+    "decay_steps":	5000,
+    "decay_rate":	0.95,
+
+    "start_pref_e":	0.02,
+    "limit_pref_e":	1,
+    "start_pref_f":	1000,
+    "limit_pref_f":	1,
+    "start_pref_v":	0,
+    "limit_pref_v":	0,
+
+    "seed":		1,
+
+    "_comment": " display and restart",
+    "_comment": " frequencies counted in batch",
+    "disp_file":	"lcurve.out",
+    "disp_freq":	100,
+    "numb_test":	50,
+    "save_freq":	100,
+    "save_ckpt":	"model.ckpt",
+    "load_ckpt":	"model.ckpt",
+    "disp_training":	true,
+    "time_training":	true,
+    "profiling":	false,
+    "profiling_file":	"timeline.json",
+
+    "_comment":		"that's all"
+}
+