Repository: deepmodeling/rid-kit (language: python)

Issue:
First refactor of Rid-kit
This is the first pr of rid-kit, the python package is named by 'rid'.
please review the structure and README.md.

unittest will be submitted in the next pr.


Discussion / hints:

Fix patch:
diff --git a/.gitignore b/.gitignore
index 8bb4cf74..a3530e79 100644
--- a/.gitignore
+++ b/.gitignore
@@ -7,4 +7,5 @@ dist
 build
 tmp*
 _version.py
-.coverage
\ No newline at end of file
+.coverage
+.DS_Store
\ No newline at end of file
diff --git a/README.md b/README.md
index 3da5e824..9f64cdac 100644
--- a/README.md
+++ b/README.md
@@ -8,6 +8,7 @@
 - [Installation of DeepMD potential](#installation-of-dp)
 - [Configure simulations](#configure-simulations)
 - [Configure machine resources](#configure-machine-resources)
+- [Configure MCMC dimension reduction](#configure-mcmc)
 - [Troubleshooting](#troubleshooting)
 
 # About Rid-kit
@@ -33,41 +34,121 @@ export BOHRIUM_PROJECT_ID="<bohrium-project-id>"
 
 ## Install Rid-kit
 Install the latest rid-kit
-```
+
+```bash
 git clone git@github.com:PKUfjh/rid-kit.git
 cd rid-kit
 git checkout dflow
 pip install setuptools_scm
 pip install .
 ```
+
 ## Run an example
 Change to the rid-kit directory
-```
+
+```bash
 cd rid-kit
 ```
+
 Run a example of Ala-dipeptide using dihedral as CVs (change to your own Bohrium account information)
-```
+
+```bash
 rid submit -i ./tests/data/000 -c ./rid/template/rid_gmx_dih.json -m ./rid/template/machine_bohrium_k8s.json
 ```
+
 You can also run the example on a Slurm machine (But you need to configure a conda enviroment on the slurm, see below)
-```
+
+```bash
 rid submit -i ./tests/data/000 -c ./rid/template/rid_gmx_dih.json -m ./rid/template/machine_slurm_k8s.json
 ```
+
 You can also run the example using distance as CVs (This will use constrained md as the mean force calculator)
-```
+
+```bash
 rid submit -i ./tests/data/000 -c ./rid/template/rid_gmx_dis.json -m ./rid/template/machine_bohrium_k8s.json
 ```
 
+You can specify the workflow name by providing WORKFLOW_ID after "-d", for example:
+
+```bash
+rid submit -i ./tests/data/000 -c ./rid/template/rid_gmx_dih.json -m ./rid/template/machine_bohrium_k8s.json -d ala-dipeptide-1
+```
+
+**Note that the defined workflow-id should only contain lower case alphanumeric character, and specifal character "-".**
+
 Note that if you want to use constrained MD as the mean force calculator, apart from setting `method` to be `constrained` in the `label_config`, you should add `[ constraints ]` line corresponding to the `[ moleculartype ]` in your input `topology` file yourself, since gromacs specifies constraints information for each `[ moleculartype ]`.
 
 Also note that, since gromacs only supports constrained MD for `distance CV`, the constrained MD simulation in rid-kit only supports `distance CV` at this moment.
 
-# Run the Workflow without k8s environment
+## Continue from an old workflow
+
+Using `resubmit` to continue from an old workflow
+
+```bash
+# suppose the original workflow id is OLD_ID
+rid resubmit -i ./tests/data/000 -c ./rid/template/rid_gmx_dis.json -m ./rid/template/machine_bohrium_k8s.json OLD_ID -d NEW_ID 
+```
+
+If you want to resubmit from a particular `iteration` and `step`:
+
+```bash
+rid resubmit -i your_dir -c path_to_rid.json -m path_to_machine.json OLD_ID -t ITERATION-ID -p STEP-KEY -d NEW_ID
+```
+
+ The `ITERATION-ID` is just `n`th iteration the workflow has been executed. The `STEP-KEY` in rid includes the following steps: `prep-exploration`, `run-exploration`, `prep-select`, `run-select`, `prep-label`, `run-label`, `collect`, `merge`, `train`, `model-devi`.
+
+## Download files from the workflow
+
+```bash
+rid download WORKFLOW_ID -p STEP-KEY -f FILE_NAME -a ITERATION_START -e ITERATION_END -o OUTPUT_DIR
+```
+
+typically we want the trajectories information from each exploration step, suppose we run a workflow for 20 iterations.
+
+```bash
+rid download WORKFLOW_ID -p run-exploration -f trajectory -a 1 -e 20 -o my_protein_out
+```
+
+`ITERATION_START`(default 1) specifies from which iteration to start download, `ITERATION_END`(default 100) specifies to which iteration to end download. `OUTPUT_DIR`(default "./") is the output directories of the download.
+
+The `STEP-KEY` in rid includes the following steps: `prep-exploration`, `run-exploration`, `prep-select`, `run-select`, `prep-label`, `run-label`, `collect`, `merge`, `train`, `model-devi`.
+
+The `FILE_NAME` in `prep-exploration` includes `task_path`.
+
+The `FILE_NAME` in `run-exploration` includes `trajectory`, `conf_out`,`plm_out`,`md_log`,`bias_fig`, `model_devi_fig`, `dp_model_devi_fig`, `dp_model_devi`, `dp_selected_indices`, `dp_selected_confs`, `projected_fig`.
+
+The `FILE_NAME` in `prep-select` includes `cluster_fig`, `cluster_selection_index`, `cluster_selection_data`.
+
+The `FILE_NAME` in `run-select` includes `selected_confs`, `selected_cv_init`, `model_devi`, `selected_indices`,`selected_conf_tags`.
+
+The `FILE_NAME` in `prep-label` includes `task_path`.
+
+The `FILE_NAME` in `run-label` includes `plm_out`, `cv_forces`, `mf_info`, `mf_fig`,`md_log`.
+
+The `FILE_NAME` in `collect` includes `data_new`.
+
+The `FILE_NAME` in `merge` includes `data_raw`.
+
+The `FILE_NAME` in `train` includes `model`, `train_log`, `train_fig`.
+
+The `FILE_NAME` in `model_devi` includes `model_devi`, `model_devi_png`.
+
+## Reduce the dimension of free energy model
+After the Rid-kit Run, the workflow will generate several numbers of free energy models (.pb). The Rid-kit currently support MCMC to reduce the dimension of the free energy model, for example:
+```bash
+rid redim -i ./test/data/models -c ./rid/template/rid_mcmc_cv_dih.json -m ./rid/template/machine_bohrium_k8s_mcmc.json
+```
+Then you will get the projected free energy surface for ala-dipeptide
+- ![image](docs/pics/mcmc_ala_1cv.png)
+- ![image](docs/pics/mcmc_ala_2cv.png)
+
+# Run the Workflow without k8s environment`
 To run the workflow without k8s environment, one can use the `Debug` mode of `Dflow`. In this mode however, one can not monitor the workflow in the `Argo` UI.
 ## Run an example
-If one wants to run the workflow on the `Slurm` machine locally, type
-```
-DFLOW_DEBUG=1 rid submit -i /tests/data/000 -c /rid/template/rid_gmx_dih.json -m /rid/template/machine_slurm_local.json
+If one wants to run the workflow on the `Slurm` machine locally, change to the rid-kit directory and type (change to your slurm configuration)
+
+```bash
+DFLOW_DEBUG=1 rid submit -i ./tests/data/000 -c ./rid/template/rid_gmx_dih.json -m ./rid/template/machine_slurm_local.json -d ala-dipeptide-1
 ```
 
 # Main procedure of RiD
@@ -92,9 +173,9 @@ A fully connected NN will be trained via sampling data. This network will genera
 
 A more detailed description of RiD is published now, please see:
 
->  [1]  Zhang, L., Wang, H., E, W.. Reinforced dynamics for enhanced sampling in large atomic and molecular systems[J]. The Journal of chemical physics, 2018, 148(12): 124113.
+> [1]  Zhang, L., Wang, H., E, W.. Reinforced dynamics for enhanced sampling in large atomic and molecular systems[J]. The Journal of chemical physics, 2018, 148(12): 124113.
 >  
->  [2]  Wang, D., Wang, Y., Chang, J., Zhang, L., & Wang, H. (2022). Efficient sampling of high-dimensional free energy landscapes using adaptive reinforced dynamics. Nature Computational Science, 2(1), 20-29.
+> [2]  Wang, D., Wang, Y., Chang, J., Zhang, L., & Wang, H. (2022). Efficient sampling of high-dimensional free energy landscapes using adaptive reinforced dynamics. Nature Computational Science, 2(1), 20-29.
 
 # Use Rid-kit
 
@@ -117,6 +198,9 @@ Installation of the DeepMD potential support can be found in
 # Configure machine resources
 - [configure machines](docs/source/rid_machine.md)
 
+# Configure mcmc 
+- [configure mcmc dimension reduction](docs/source/mcmc_configuration.md)
+
 # Workflow Synopsis
 
 - ![image](docs/pics/rid_workflow.jpg)
diff --git a/docs/pics/mcmc_ala_1cv.png b/docs/pics/mcmc_ala_1cv.png
new file mode 100644
index 00000000..437bfa44
Binary files /dev/null and b/docs/pics/mcmc_ala_1cv.png differ
diff --git a/docs/pics/mcmc_ala_2cv.png b/docs/pics/mcmc_ala_2cv.png
new file mode 100644
index 00000000..88df327c
Binary files /dev/null and b/docs/pics/mcmc_ala_2cv.png differ
diff --git a/docs/requirements.txt b/docs/requirements.txt
index ea1ac232..4701a2d6 100755
--- a/docs/requirements.txt
+++ b/docs/requirements.txt
@@ -2,10 +2,13 @@ pydflow
 numpy
 scikit-learn
 mdtraj
-tensorflow
+dpdata
+matplotlib
+# tensorflow
 ipython
 jinja2==3.0.1
 docutils==0.17.1
 myst-parser
 m2r2
-nbsphinx
\ No newline at end of file
+nbsphinx
+sphinx_rtd_theme
\ No newline at end of file
diff --git a/docs/source/conf.py b/docs/source/conf.py
index 47962629..d97cb311 100755
--- a/docs/source/conf.py
+++ b/docs/source/conf.py
@@ -25,6 +25,8 @@
 # The full version, including alpha/beta/rc tags
 release = '2.0'
 
+nbsphinx_allow_errors = True
+
 
 # -- General configuration ---------------------------------------------------
 
diff --git a/docs/source/index.rst b/docs/source/index.rst
index 4ba2409e..b57a0aca 100755
--- a/docs/source/index.rst
+++ b/docs/source/index.rst
@@ -25,6 +25,7 @@ Contents
    Installation with DP potential <install_dp>
    Configure simulations <rid_configuration>
    Configure machine resources <rid_machine>
+   Configure mcmc dimension reduction <mcmc_configuration>
    Rid API <api/api>
    Trouble Shooting <troubleshooting/index> 
 
diff --git a/docs/source/install.md b/docs/source/install.md
index 40c897ed..5de40126 100644
--- a/docs/source/install.md
+++ b/docs/source/install.md
@@ -51,15 +51,27 @@ export PLUMED_KERNEL=${CONDA_PREFIX}/lib/libplumedKernel.@SOEXT@
 ### Install gromacs 2022.4
 Note that in order to compile gromacs with GPU support, you need to have CUDA toolkit installed. In slurm enviroment, if the software is managed by `Enviroment Module`, one can load CUDA toolkits such as `module load cuda/11.4`. Also note that in most slurm environment, there is no network connection at computing node, but if you want to compile gromacs without fftw installed, you will have to compile it on login-in node where network connection is on.
 
+However if the network is not on in the slurm machine, you can first compile fftw manually by
+```bash
+cd fftw-3.3.8
+./configure --prefix=$CONDA_PREFIX
+make
+make install
+```
+
 ```bash
 tar -xzvf gromacs-2022.4.tar.gz
 cd gromacs-2022.4
 plumed patch -p # make sure you use the correct gromacs version
 mkdir build
 cd build
+# if the network is on
 cmake .. -DGMX_BUILD_OWN_FFTW=ON \
                 -DCMAKE_INSTALL_PREFIX=$CONDA_PREFIX \
                 -DGMX_GPU=CUDA
+# if the network is off
+cmake .. -DCMAKE_INSTALL_PREFIX=$CONDA_PREFIX \
+            -DGMX_GPU=CUDA
 make -j 4
 make install
 ```
diff --git a/docs/source/mcmc_configuration.md b/docs/source/mcmc_configuration.md
new file mode 100644
index 00000000..1015c2b6
--- /dev/null
+++ b/docs/source/mcmc_configuration.md
@@ -0,0 +1,85 @@
+# Setting Parameters of MCMC dimension reduction
+
+`RiD-kit` uses a JSON-format file to configure MCMC dimension reduction. Here we explain these parameters one by one.
+
+```JSON
+{
+        "name": "rid-mcmc",
+        "init_models": ["model_000.pb","model_001.pb","model_002.pb","model_003.pb"],
+        
+        "MCMC_Config": {
+            "cv_dimension": 2,
+            "numb_steps": 10000,
+            "numb_walkers": 500,
+            "temperature": 300,
+            "proj_info": {
+                "proj_mode": "cv",
+                "proj_cv_index": [0,1]
+            },
+            "cv_type": "dih",
+            "bins": 101
+        }
+    
+},
+{
+        "name": "rid-mcmc",
+        "init_models": ["model_000.pb","model_001.pb","model_002.pb","model_003.pb"],
+        
+        "MCMC_Config": {
+            "cv_dimension": 20,
+            "numb_steps": 100000,
+            "numb_walkers": 2000,
+            "cv_upper_bound": [4.5,5.0,6.0,4.2,4.0,4.0,3.5,4.0,4.0,3.5,3.2,3.0,3.5,3.5,4.0,4.5,4.5,3.75,3.8,4.0],
+            "cv_lower_bound": [1.0,2.0,3.5,2.8,2.5,2.0,1.2,1.5,1.5,1.0,1.8,0.5,1.0,1.0,1.5,1.5,1.5,1.75,2.2,2.6],
+            "proj_info": {
+                "proj_mode": "cv",
+                "proj_cv_index": [0,1]
+            },
+            "cv_type": "dis",
+            "bins": 101
+        }
+    
+}
+```
+* **`"init_models"`** `(List[str])` List of trained free energy models generated in the RiD run.
+* **`"cv_dimension"`** `int` An integer to represent the number of CVs used in RiD.
+* **`"numb_steps"`** `int` An integer to represent the number of steps for MCMC run.
+* **`"numb_walkers"`** `int` Number of parallel walkers used in the MCMC run.
+* **`"temperature"`**  `int` Temperature used in the RiD run and MCMC simulation.
+* **`"cv_upper_bound"`** `Optional[(List[float])]` Upper bound for each dimension of CVs, this is usually set for distance CVs.
+* **`"cv_lower_bound"`** `Optional[(List[float])]` Lower bound for each dimension of CVs, this is usually set for distance CVs.
+* **`"proj_mode"`** `str` The mode for projecting 2D CV, support `cv` and `path` mode. `cv` mode representing the common projection on the selected CV index, `path` mode representing the projection on the user defined path CV.
+* **`"proj_cv_index"`** `(List[int])` Projected index for the 2D CV. Note that 1D projection is done for all CVs.
+* **`"cv_type"`** `str` CV type for projection, support `dih` and `dis`, representing dihedral CVs and distance CVs.
+* **`"bins"`** `int` Bins for projecting CVs at each dimension.
+
+Another example is given to project 2D path CV.
+```JSON
+{
+        "name": "rid-mcmc",
+        "init_models": ["model_000.pb","model_001.pb","model_002.pb","model_003.pb"],
+        
+        "MCMC_Config": {
+            "cv_dimension": 13,
+            "numb_steps": 100000,
+            "numb_walkers": 2000,
+            "cv_upper_bound": [5.33507,4.43402,4.59772,3.40287,4.67268,5.12813,3.28143,3.74967,3.58144,5.34143,7.13014,4.58024,7.41678],
+            "cv_lower_bound": [2.38722,0.99005,3.99886,2.08580,3.81799,3.93331,2.768175,3.145618,3.161886,4.661101,5.286805,4.076754,2.893047],
+            "proj_info": {
+                "proj_mode": "path",
+                "proj_cv_index": [0,1,2,3],
+                "path_list": [[2.38722,0.99005,3.99886,2.08580],[2.64067,1.20917,3.73831,1.86004],[2.90639,1.50430,3.47935,1.65185],
+                [3.18130,1.83919,3.22240,1.46871],[3.46322,2.19574,2.96795,1.32108],[3.75057,2.56491,2.71672,1.22191],
+                [4.04218,2.94198,2.46968,1.18344],[4.33720,3.32424,2.22824,1.21147],[4.63498,3.71010,1.99443,1.30170]],
+                "path_lm": 0.25
+            },
+            "cv_type": "dis",
+            "bins": 101
+        }
+    
+}
+```
+In the path CV projection mode, the additonal parameters are:
+* **`"proj_cv_index"`** `List[int]` Set of CV index used to defined the path CV.
+* **`"path_list"`** `List[float]` List of CVs to defined the path.
+* **`"path_lm"`** `float` Lambda in the path CV definition.
\ No newline at end of file
diff --git a/docs/source/troubleshooting/installation.md b/docs/source/troubleshooting/installation.md
index becad95b..e08cee62 100644
--- a/docs/source/troubleshooting/installation.md
+++ b/docs/source/troubleshooting/installation.md
@@ -10,6 +10,16 @@ add-apt-repository -y ppa:ubuntu-toolchain-r/test
 apt-get update
 apt-get install -y libstdc++6
 ```
+Or in slurm machine you might encounter this error in conda enviroment
+```
+/data1/ddwang/software/anaconda3/envs/rid-drop/lib/libgrpc.so.20: undefined reference to `std::__throw_bad_array_new_length()@GLIBCXX_3.4.29'
+/data1/ddwang/software/anaconda3/envs/rid-drop/lib/libtensorflow_framework.so.2: undefined reference to `std::__cxx11::basic_ostringstream<char, std::char_traits<char>, std::allocator<char> >::basic_ostringstream()@GLIBCXX_3.4.26'
+/data1/ddwang/software/anaconda3/envs/rid-drop/lib/libtensorflow_cc.so.2: undefined reference to `std::__cxx11::basic_stringstream<char, std::char_traits<char>, std::allocator<char> >::basic_stringstream()@GLIBCXX_3.4.26'
+```
+Try to install libgcc in conda
+```
+conda install libgcc
+```
 
 **You may encounter the error when compiling gromacs**
 ```
diff --git a/docs/source/tutorial.ipynb b/docs/source/tutorial.ipynb
index 8e164da9..161d52fb 100644
--- a/docs/source/tutorial.ipynb
+++ b/docs/source/tutorial.ipynb
@@ -264,7 +264,39 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "! rid submit -i /tests/data/000 -c /rid/template/rid_gmx_dih.json -m /rid/template/machine_bohrium_local.json"
+    "! rid submit -i /tests/data/000 -c /rid/template/rid_gmx_dih.json -m /rid/template/machine_bohrium_k8s.json"
+   ]
+  },
+  {
+   "attachments": {},
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "You can specify the workflow name by providing WORKFLOW_ID after \"-d\", for example:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "! rid submit -i ./tests/data/000 -c ./rid/template/rid_gmx_dih.json -m ./rid/template/machine_bohrium_k8s.json -d ala-dipeptide-1"
+   ]
+  },
+  {
+   "attachments": {},
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "**Note that the defined workflow-id should only contain lower case alphanumeric character, and specifal character \"-\".**"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Note that the defined workflow-id should only contain lower case alphanumeric character, and specifal character \"-\"."
    ]
   },
   {
@@ -331,7 +363,8 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "! rid resubmit -i your_dir -c path_to_rid.json -m path_to_machine.json Workflow-ID"
+    "# suppose the original workflow id is OLD_ID\n",
+    "! rid resubmit -i your_dir -c path_to_rid.json -m path_to_machine.json OLD_ID -d NEW_ID"
    ]
   },
   {
@@ -347,7 +380,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "! rid resubmit -i your_dir -c path_to_rid.json -m path_to_machine.json Workflow-ID -t ITERATION-ID -p STEP-KEY"
+    "! rid resubmit -i your_dir -c path_to_rid.json -m path_to_machine.json OLD_ID -t ITERATION-ID -p STEP-KEY -d NEW_ID"
    ]
   },
   {
diff --git a/docs/source/tutorial.md b/docs/source/tutorial.md
new file mode 100644
index 00000000..264c1150
--- /dev/null
+++ b/docs/source/tutorial.md
@@ -0,0 +1,195 @@
+# Tutorial
+
+In this tutorial, we will learn how to deploy `RiD-kit` with `dflow` and `kubenete` and run a simple case of alanine dipeptide.
+
+
+## Installation of `dflow` and `rid-kit`
+With the power of `dflow`, users can easily minitor the whole workflow of RiD tasks and dispatch their tasks to various computational resources. Before you use it, you should have `dflow` installed on your host computer (your PC or a remote server). 
+
+It it necessary to emphasize that, the computational nodes and monitor nodes are seperated. With `dflow`, you can deploy `dflow` and `rid` on your PC and achieve expensive computation on other resources (like `Slurm` and Cloud Platform) without any further effort.
+
+Instructions of `dflow` installation are provided in detail on its [Github page](https://github.com/deepmodeling/dflow#Installdflow). Prerequisites of `dflow` usage are `Docker` and `Kubenetes`, where their main pages ([Docker](https://docs.docker.com/engine/install/) & [Kubenetes](https://kubernetes.io/docs/tasks/tools/)) include how you can install them. Besides, `dflow` repo also provides with easy-install shell scripts on [dflow/scripts](https://github.com/deepmodeling/dflow/tree/master/scripts) to install `Docker` & `Kubenetes` & `dflow` and make port-forwarding.
+
+Instead of using the easy-install scripts, it is recommended to configure your own minikube enviroment
+
+Installation of `minikube` is very easy, just follow the page [Minikube](https://minikube.sigs.k8s.io/docs/start/) based on your machine archetecture.
+
+Then use minikube to start a k8s environment. This can be done by just type `minikube start`, but since running rid-kit requires relatively huge amount of resources, it is recommended to allocate more resources when executing `minikube start`, just as follows.
+
+> **Note:**
+>
+> Don't try to run `minikube` with root privileges, otherwise an error may occur:
+>
+> `Exiting due to DRV_AS_ROOT: The "docker" driver should not be used with root privileges.`
+>
+
+```bash
+#choose the location of the minikube enviroment
+export MINIKUBE_HOME=~/.minikube  
+#allocate enough memory based on your machine in case of high parallelism
+minikube start --cpus 8 --memory 8192mb --kubernetes-version=1.23.9 --image-mirror-country='cn'
+# mount the storage path if you are using machine with shared memory storage system, change the minio host path accordingly
+minikube start --cpus 8 --memory 8192mb --kubernetes-version=1.23.9 --mount --mount-string="/path_on_your_machine:/data2" --image-mirror-country='cn'
+```
+
+A further step to  configure `argo` service is to run:
+
+```bash
+kubectl create ns argo
+kubectl apply -n argo -f https://raw.githubusercontent.com/deepmodeling/dflow/master/manifests/quick-start-postgres.yaml
+```
+
+If you want to mount your local path into the minio path, change the minio host path in the yaml file
+
+```yaml
+      - hostPath:
+          path: /data2/minio
+          type: DirectoryOrCreate
+```
+
+Now you should have installed `Docker` and `minikube` properly. Run commands to check their status. For `minikube`, you should wait util all servers keep `running`. This may take a couple of minutes.
+
+```bash
+# check minikube status
+minikube status
+# check k8s status in argo namespace
+kubectl get pods -n argo
+```
+
+### Installation of RiD-kit
+
+Now we install `rid-kit` on the host machine. To meet the minimum requirments, the prerequisites of third-party python package should be installed:
+
+* tensorflow-cpu or gpu
+* mdtraj
+* numpy
+* scikit-learn
+* pydflow
+* google
+* mdtraj
+* dpdata
+* parmed
+* dpdispatcher
+* lbg
+* matplotlib
+
+which are also listed in `rid-kit/requirements.txt`. Then change directory to `rid-kit` repo and run:
+
+```bash
+# the rid-kit repo path
+git checkout dflow
+pip install setuptools_scm
+pip install .
+```
+
+## Configuration of Computational Environment
+
+In RiD workflow, `dflow` helps send computation tasks to resources with peoper environment configured. 
+
+There are four main modules and several workflow steps in RiD procedures and each module or step needs different environments:
+
+* `Exploration/Sampling`:  `Gromacs`, `PLUMED2` modified by `DeepFE.cpp`, `Tensorflow` C++ interface. (prefer GPU)
+* `Selection`:  `Tensroflow` Python interface.
+* `Labeling`:  `Gromacs`, `PLUMED2`. (prefer GPU)
+* `Training`:  `Tensroflow` Python interface. (prefer GPU).
+* `Workflow steps`:  `Python`.
+
+`dflow` supports different resources including `Slurm` clusters, `K8S` local machines and `Cloud Server`.
+
+* For `Slurm`, configure computational environments on your `Slurm` following the instructions in [Environment settings](https://github.com/PKUfjh/rid-kit/blob/dflow/docs/source/install.md). With `dflow`, `rid-kit` send tasks to `Slurm` nodes from the host machines remotely without manually logging in the cluster.
+* For local resources, just use the docker images we have built. No further manual configuration needed. We also provide `Dockerfile` of our images to enable flexible modification.
+* For `Cloud Server`, like `Bohrium`, use public images and no further manual configuration needed.
+
+We highly recommend using `Bohrium` to do computation since no enviroment installation is required, you only need to use the public image.Specifically, the public images for rid-kit involves the following: `registry.dp.tech/public/pkufjhdocker/rid-gmx-exploration:stable` for doing `Exploration`, `registry.dp.tech/public/pkufjhdocker/rid-gmx-plumed:stable` for doing `Labeling`, `registry.dp.tech/public/pkufjhdocker/rid-gmx-tf:stable` for doing `Selection`, `registry.dp.tech/public/pkufjhdocker/rid-tf-gpu:stable` for `Training` and `registry.dp.tech/public/pkufjhdocker/rid-tf-cpu:stable` for doing other steps.
+
+## Prepare machine configuration JSON.
+
+`rid-kit` uses `JSON` file to manage resources. In `machine.json`, define your own `resources` and dispatch `tasks` to them.
+* In key `resources`, you define your own resources types. Resource names and their numbers are custom.
+* In key `tasks`, you distribute resources you have defined to tasks of RiD. Do not change task names in `tasks` as they are fixed in codes.
+
+Generally, we would like to run low-cost tasks on cpu nodes or locally and submit high-cost tasks to `Slurm` or `Clouds`. 
+
+If submitting to mixed `Slurm` and `local` environment, a `machine.json` example is in [rid/template/machine_slurm_k8s.json](https://github.com/PKUfjh/rid-kit/blob/dflow/rid/template/machine_slurm_k8s.json):
+
+If submitting to `Bohrium` cloud platform, a `machine.json` example is in [rid/template/machine_bohrium.json](https://github.com/PKUfjh/rid-kit/blob/dflow/rid/template/machine_bohrium.json).
+
+If submitting to mixed `Bohrium` and `local` platform, a `machine.json` example is in [rid/template/machine_bohrium_k8s.json](https://github.com/PKUfjh/rid-kit/blob/dflow/rid/template/machine_bohrium_k8s.json).
+
+If submitting to `Slurm` environment without `k8s` enviroment, a `machine.json` example is in [rid/template/machine_slurm_local.json](https://github.com/PKUfjh/rid-kit/blob/dflow/rid/template/machine_slurm_local.json):
+
+If you submit jobs to `Slurm` enviroment, you will have to compile the computation enviroment on `Slurm` server for yourself. While in `Bohrium` cloud platform the computation images are all set, no extra compilation is needed. So we recommend to submit jobs to `Cloud` enviroment like `Bohrium`.
+
+For more details about setting the machine resources, check [Rid machine resourse](https://github.com/PKUfjh/rid-kit/blob/dflow/docs/source/rid_machine.md)
+
+## Prepare rid configuration JSON.
+If you want to use `gromacs` as sampler, [rid/template/rid_gmx_dih.json](https://github.com/PKUfjh/rid-kit/blob/dflow/rid/template/rid_gmx_dih.json) is an example.
+
+If you want to use `gromacs` with `DeepMD` potential as sampler, [rid/template/rid_gmx_dp.json](https://github.com/PKUfjh/rid-kit/blob/dflow/rid/template/rid_gmx_dp.json) is an example.
+
+If you want to use `lammps` with `DeepMD` potential as sampler, [rid/template/rid_lmp_dp.json](https://github.com/PKUfjh/rid-kit/blob/dflow/rid/template/rid_lmp_dp.json) is an example.
+
+For more details about setting the rid configuration, check [Rid configuration](https://github.com/PKUfjh/rid-kit/blob/dflow/docs/source/rid_configuration.md)
+
+## Get Started!
+
+Assume you have learn the basic knowledge of reinforced dynamics which we won't describe again here. 
+
+Users can monitor workflows from browser UI. To enable that, you should forward ports of `argo` and `minio`. These could be achieved by `rid port-forward`.
+
+```bash
+rid port-forward
+```
+
+In this case, we try to explore the phase space of alanine dipeptide. Prepare your initial conformation files in `.gro` format, topology file in `.top` format and configuration file `rid.json`. For convenience, we have prepared on at `rid/template/rid.json`. 
+Remember also provide your own forcefield files. Collect all these files into a directory and feed its path to `rid-kit` by flag `-i`.
+
+A minimum case was prepared in `rid-kit/tests/data/000`. After configurating your `machine.json`, then run `rid submit`:
+
+```bash
+rid submit -i ./tests/data/000 -c ./rid/template/rid_gmx_dih.json -m ./rid/template/machine_bohrium_k8s.json
+```
+
+You can specify the workflow name by providing WORKFLOW_ID after "-d", for example:
+
+```bash
+rid submit -i ./tests/data/000 -c ./rid/template/rid_gmx_dih.json -m ./rid/template/machine_bohrium_k8s.json -d ala-dipeptide-1
+```
+
+**Note that the defined workflow-id should only contain lower case alphanumeric character, and specifal character "-".**
+
+INFO indicates that this task has been submitted succussfully. Record this workflow ID as we may use it later.
+
+Visit the `url` given by the last two lines, all workflows and corresponding files are listed on UI.
+
+Command lines are also supported. Run `rid ls` to list your workflows and their status.
+
+```bash
+! rid ls
+```
+
+`rid-kit` is based on `dflow`, `argo` and `minikube`. So further complex and flexible managements of workflows can be achieved by their command lines. like `kubectl get pods -n argo` and `argo show`.
+
+For failed tasks, you may want to remove them or resubmit them from the failure steps.
+
+For `remove`:
+
+```bash
+# rid rm task-ID
+rid rm reinforced-dynamics-bsc7j 
+```
+
+For `resubmit` to modify and continue workflow:
+
+```bash
+# suppose the original workflow id is OLD_ID
+rid resubmit -i your_dir -c path_to_rid.json -m path_to_machine.json OLD_ID -d NEW_ID
+```
+
+If you want to resubmit from a particular `iteration` and `step`:
+
+```bash
+rid resubmit -i your_dir -c path_to_rid.json -m path_to_machine.json OLD_ID -t ITERATION-ID -p STEP-KEY -d NEW_ID
+```
+
+The `Workflow-ID` is something like `reinforced-dynamics-jq4jn` appeared in the argo UI. If the workflow is archived, its name will appear as somethinig like `a8463748-e15a-4f2c-882b-bfd981a76dac`. If you want to resubmit a archived workflow, you have to provide the archived ID rather than the initial ID. The `ITERATION-ID` is just `n`th iteration the workflow has been executed. The `STEP-KEY` in rid includes the following steps: `prep-exploration`, `run-exploration`, `prep-select`, `run-select`, `prep-label`, `run-label`, `post-label`, `collect`, `merge`, `train`.
diff --git a/install/lmp_exploration_dp.dockerfile b/install/lmp_exploration_dp.dockerfile
index 70393ca2..d51892a2 100644
--- a/install/lmp_exploration_dp.dockerfile
+++ b/install/lmp_exploration_dp.dockerfile
@@ -30,7 +30,7 @@ RUN conda init bash && source /root/.bashrc \
     && conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/pkgs/main/  \
     && conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/pytorch/  \
     && conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/conda-forge  \
-    && conda install python=3.9 libtensorflow_cc=2.6.2=*cuda110* dpdata nccl numpy -y
+    && conda install python=3.9 libtensorflow_cc=2.6.2=*cuda110* tensorflow=2.6.2=*cuda110* dpdata nccl numpy -y
 
 # Solve library inconsistency
 RUN source activate base && rm ${CONDA_PREFIX}/lib/libtinfo.so* && ln -s /usr/lib/x86_64-linux-gnu/libtinfo.so.6 ${CONDA_PREFIX}/lib/libtinfo.so.6 \
diff --git a/rid/__init__.py b/rid/__init__.py
index acc3d8ba..405deae5 100644
--- a/rid/__init__.py
+++ b/rid/__init__.py
@@ -11,4 +11,5 @@
 SRC_ROOT = __path__[0]
 config["extender_image_pull_policy"] = "IfNotPresent"
 config["util_image_pull_policy"] = "IfNotPresent"
-config["dispatcher_image_pull_policy"] = "IfNotPresent"
\ No newline at end of file
+config["dispatcher_image_pull_policy"] = "IfNotPresent"
+config["save_keys_in_global_outputs"] = False
\ No newline at end of file
diff --git a/rid/common/mol.py b/rid/common/mol.py
index bc9f070c..c8152ad8 100644
--- a/rid/common/mol.py
+++ b/rid/common/mol.py
@@ -3,7 +3,6 @@
 from typing import List, Dict, Sequence, Union
 import logging
 import mdtraj as md
-import dpdata
 from mdtraj.geometry.dihedral import _atom_sequence, PHI_ATOMS, PSI_ATOMS
 import numpy as np
 from rid.constants import sel_gro_name_gmx, sel_gro_name
@@ -122,32 +121,5 @@ def slice_xtc(
             selected_idx = selected_idx,
             output_format=output
         )
-    else:
-        raise RuntimeError("Unknown Style for Slicing Trajectory.")
-
-def slice_dump_dpdata(
-    dump: str,
-    walker_idx,
-    selected_idx,
-    output_format:str
-):
-    traj = dpdata.System(dump, fmt="lammps/dump")
-    for sel in selected_idx:
-        traj.to('lammps/lmp', output_format.format(walker=walker_idx,idx=sel), frame_idx=sel)
-
-def slice_dump(
-    dump: str,
-    walker_idx: int,
-    selected_idx,
-    output:str,
-    style: str = "dpdata"
-):
-    if style == "dpdata":
-        slice_dump_dpdata(
-            dump = dump,
-            walker_idx= walker_idx,
-            selected_idx = selected_idx,
-            output_format=output
-        )
     else:
         raise RuntimeError("Unknown Style for Slicing Trajectory.")
\ No newline at end of file
diff --git a/rid/common/mol_dpdata.py b/rid/common/mol_dpdata.py
new file mode 100644
index 00000000..393315de
--- /dev/null
+++ b/rid/common/mol_dpdata.py
@@ -0,0 +1,41 @@
+import dpdata
+from typing import List, Dict, Sequence, Union
+
+def slice_dump_dpdata(
+    dump: str,
+    walker_idx,
+    selected_idx,
+    output_format:str,
+    type_map: List = [],
+):
+    if type_map != []:
+        traj = dpdata.System(dump, fmt="lammps/dump", type_map=type_map)
+    else:
+        traj = dpdata.System(dump, fmt="lammps/dump")
+    if output_format.endswith("lmp"):
+        for sel in selected_idx:
+            traj.to('lammps/lmp', output_format.format(walker=walker_idx,idx=sel), frame_idx=sel)
+    elif output_format.endswith("gro"):
+        for sel in selected_idx:
+            traj.to('gromacs/gro', output_format.format(walker=walker_idx,idx=sel), frame_idx=sel)
+    else:
+        raise ValueError("Invalid output format, only support gmx and lmp")
+
+def slice_dump(
+    dump: str,
+    walker_idx: int,
+    selected_idx,
+    output:str,
+    type_map: List = [],
+    style: str = "dpdata"
+):
+    if style == "dpdata":
+        slice_dump_dpdata(
+            dump = dump,
+            walker_idx= walker_idx,
+            selected_idx = selected_idx,
+            type_map = type_map,
+            output_format=output
+        )
+    else:
+        raise RuntimeError("Unknown Style for Slicing Trajectory.")
\ No newline at end of file
diff --git a/rid/constants.py b/rid/constants.py
index cf658bec..2f51c377 100644
--- a/rid/constants.py
+++ b/rid/constants.py
@@ -1,13 +1,14 @@
 # rid task name format
 argo_namespace = "argo"
 walker_tag_fmt = "{idx:03d}"
-init_conf_name = "conf_{idx:03d}.gro"
+init_conf_gmx_name = "conf_{idx:03d}.gro"
+init_conf_lmp_name = "conf_{idx:03d}.lmp"
 init_input_name = "input_{idx:03d}.lammps"
 explore_task_pattern = "{:03d}"
 explore_task_file = "explore_{walker:03d}.pkl"
 cluster_selection_data_name = "cls_sel.out.npy"
 cluster_selection_index_name = "cls_sel.ndx.npy"
-sel_ndx_name = "sel.ndx"
+sel_ndx_name = "sel.ndx.npy"
 cv_init_label = "cv_init_{walker:03d}_{idx:d}.out"
 model_devi_name = "model_devi.txt"
 label_task_pattern = "{:03d}"
@@ -58,6 +59,7 @@
 
 # Dp config file
 dp_config_name = "dp_config"
+dp_sel_ndx = "dp_sel.ndx"
 
 # figure file name
 mf_fig = "mf_average.png"
@@ -65,13 +67,23 @@
 model_devi_fig = "model_devi.png"
 cluster_fig = "cluster.png"
 train_fig = "train_{tag}.png"
+dp_model_devi_fig = "dp_model_devi.png"
+
+# MCMC file name
+mcmc_1cv_dir_name = "mcmc_1cv"
+mcmc_2cv_name = "mcmc_2cv.dat"
+mcmc_1cv_name = "mcmc_1cv_{tag}.dat"
+mcmc_2cv_fig = "mcmc_2cv.png"
+mcmc_1cv_fig = "mcmc_1cv_{tag}.png"
+mcmc_2cv_fig_separate = "mcmc_2cv_{tag}.png"
 
 # Units
-kb = 8.617333E-5
+kb = 8.617333E-5 # ev
 kbT = (8.617333E-5) * 300
 beta = 1.0 / kbT
-f_cvt = 96.485
+f_cvt = 96.485 # ev to kj/mol
 inverse_f_cvt = 1 / f_cvt
+kcal2kj = 4.184
 
 # precision
 model_devi_precision = "%.6e"
\ No newline at end of file
diff --git a/rid/entrypoint/download.py b/rid/entrypoint/download.py
new file mode 100644
index 00000000..680d2ff8
--- /dev/null
+++ b/rid/entrypoint/download.py
@@ -0,0 +1,45 @@
+from dflow import Workflow, download_artifact
+import numpy as np
+import os
+from datetime import datetime
+
+from dflow import config, s3_config
+from dflow.plugins import bohrium
+from dflow.plugins.bohrium import TiefblueClient
+
+
+def download_rid(workflow_id, pod, file, iteration_start, iteration_end, outputs):
+    wf = Workflow(id=workflow_id)
+    all_steps = wf.query_step()
+
+    iter_a = int(iteration_start)
+    iter_e = int(iteration_end)+1
+    print("start from iteration",iter_a)
+    print("end to iteration",iter_e-1)
+    for step in all_steps:
+        if step["type"] == "Pod":
+            for index in range(iter_a, iter_e):
+                # print(step["key"])
+                if step['key'] == "iter-%03d-%s-0"%(index,pod):
+                    print("downloading "+step["key"])
+                    download_artifact(step.outputs.artifacts[file], path="%s/iter%s"%(outputs,index))
+                elif step['key'] == "iter-%03d-%s"%(index,pod):
+                    print("downloading "+step["key"])
+                    download_artifact(step.outputs.artifacts[file], path="%s/iter%s"%(outputs,index))
+
+# with open("time_profile.txt","w") as f:
+#     for step in all_steps:
+#         if step['key'] is not None:
+#             f.write("step key:"+str(step["key"]))
+#         if step['startedAt'] is not None:
+#             start_time_str = str(step['startedAt'])
+#             f.write(" start time:"+str(start_time_str))
+#             start_time = datetime.strptime(start_time_str, "%Y-%m-%dT%H:%M:%SZ")
+#         if step['finishedAt'] is not None:
+#             end_time_str = str(step['finishedAt'])
+#             f.write(" finish time:"+str(end_time_str))
+#             end_time = datetime.strptime(end_time_str, "%Y-%m-%dT%H:%M:%SZ")
+#             # Calculate the duration of the task
+#             duration = end_time - start_time
+#             f.write(" duration:"+str(duration))
+#         f.write("\n")
\ No newline at end of file
diff --git a/rid/entrypoint/main.py b/rid/entrypoint/main.py
index fccb7cc2..e666b22d 100644
--- a/rid/entrypoint/main.py
+++ b/rid/entrypoint/main.py
@@ -29,6 +29,9 @@
 from .resubmit import resubmit_rid
 from .label import label_rid
 from .relabel import relabel_rid
+from .redim import redim_rid
+from .reredim import reredim_rid
+from .download import download_rid
 from .info import information
 from .server import forward_ports
 from .cli import rid_ls, rid_rm
@@ -107,6 +110,9 @@ def main_parser() -> argparse.ArgumentParser:
     parser_run.add_argument(
         "--machine", "-m", help="Machine configuration.", dest="machine"
     )
+    parser_run.add_argument(
+        "--id", "-d", help="Workflow user defined ID", default = None, dest="id"
+    )
 
     # resubmit
     parser_rerun = subparsers.add_parser(
@@ -126,6 +132,9 @@ def main_parser() -> argparse.ArgumentParser:
     parser_rerun.add_argument(
         "--machine", "-m", help="Machine configuration.", dest="machine"
     )
+    parser_rerun.add_argument(
+        "--id", "-d", help="Workflow user defined ID", default = None, dest="id"
+    )
     parser_rerun.add_argument(
         "--iteration", "-t", help="restart from t-th iteration.", default = None, dest="iteration"
     )
@@ -204,13 +213,65 @@ def main_parser() -> argparse.ArgumentParser:
         formatter_class=argparse.ArgumentDefaultsHelpFormatter,
     )
     parser_redim.add_argument(
-        "--networks", "-n", help="Training data."
+        "--mol", "-i", help="Neural networks path", dest="mol",
+    )
+    parser_redim.add_argument(
+        "--config", "-c", help="RiD configuration.", dest="config"
     )
     parser_redim.add_argument(
-        "--dim1", help="Dimension 1."
+        "--machine", "-m", help="Machine configuration.", dest="machine"
     )
     parser_redim.add_argument(
-        "--dim2", help="Dimension 2."
+        "--id", "-d", help="Workflow user defined ID", default = None, dest="id"
+    )
+    
+    # resubmit NN dimension reduction.
+    parser_reredim = subparsers.add_parser(
+        "reredim",
+        help="NN dimension reduction by Monte Carlo integral.",
+        formatter_class=argparse.ArgumentDefaultsHelpFormatter,
+    )
+    parser_reredim.add_argument(
+        "--mol", "-i", help="Neural networks path", dest="mol",
+    )
+    parser_reredim.add_argument(
+        "--config", "-c", help="RiD configuration.", dest="config"
+    )
+    parser_reredim.add_argument(
+        "--machine", "-m", help="Machine configuration.", dest="machine"
+    )
+    parser_reredim.add_argument(
+        "--id", "-d", help="Workflow user defined ID", default = None, dest="id"
+    )
+    parser_reredim.add_argument(
+        "WORKFLOW_ID", help="Workflow ID."
+    )
+    parser_reredim.add_argument(
+        "--pod", "-p", help="restart from the pod.", default = None, dest="pod"
+    )
+    # download
+    parser_download = subparsers.add_parser(
+        "download",
+        help="Download files in RiD workflow",
+        formatter_class=argparse.ArgumentDefaultsHelpFormatter,
+    )
+    parser_download.add_argument(
+        "WORKFLOW_ID", help="Workflow ID."
+    )
+    parser_download.add_argument(
+        "--pod", "-p", help="The pod in RiD workflow", dest="pod"
+    )
+    parser_download.add_argument(
+        "--file", "-f", help="The file to download inside the step", dest = "file"
+    )
+    parser_download.add_argument(
+        "--iteration_start", "-a", help="download from a-th iteration.", default = 1, dest="iteration_start"
+    )
+    parser_download.add_argument(
+        "--iteration_end", "-e", help="download to e-th iteration.", default = 100, dest="iteration_end"
+    )
+    parser_download.add_argument(
+        "--outputs", "-o", help="output directory of downloading", default = "./", dest="outputs"
     )
 
     # --version
@@ -275,8 +336,8 @@ def parse_submit(args):
 
 
 def log_ui():
-    logger.info('The task is displayed on "https://127.0.0.1:2746".')
-    logger.info('Artifacts (Files) are listed on "https://127.0.0.1:9001".')
+    if os.getenv("DFLOW_HOST") is not None:
+        logger.info('The task is displayed on %s.'%os.getenv("DFLOW_HOST"))
 
 
 def main():
@@ -289,6 +350,7 @@ def main():
             topology = top_file,
             rid_config = args.config,
             machine_config = args.machine,
+            workflow_id_defined = args.id,
             models = models,
             forcefield = forcefield,
             index_file = index_file,
@@ -306,6 +368,7 @@ def main():
             topology = top_file,
             rid_config = args.config,
             machine_config = args.machine,
+            workflow_id_defined = args.id,
             iteration = args.iteration,
             pod = args.pod,
             models = models,
@@ -350,6 +413,38 @@ def main():
             otherfiles = otherfiles
         )
         log_ui()
+    elif args.command == "redim":
+        logger.info("Preparing MCMC ...")
+        confs, top_file, models, forcefield, index_file, data_file, dp_files, otherfiles = parse_submit(args)
+        redim_rid(
+            rid_config = args.config,
+            machine_config = args.machine,
+            models = models,
+            workflow_id_defined = args.id,
+        )
+        log_ui()
+    elif args.command == "reredim":
+        logger.info("Preparing MCMC ...")
+        confs, top_file, models, forcefield, index_file, data_file, dp_files, otherfiles = parse_submit(args)
+        reredim_rid(
+            workflow_id=args.WORKFLOW_ID,
+            rid_config = args.config,
+            machine_config = args.machine,
+            workflow_id_defined = args.id,
+            models = models,
+            pod = args.pod
+        )
+        log_ui()
+    elif args.command == "download":
+        logger.info("Downloading files ...")
+        download_rid(
+            workflow_id = args.WORKFLOW_ID,
+            pod = args.pod,
+            file = args.file,
+            iteration_start = args.iteration_start,
+            iteration_end = args.iteration_end,
+            outputs = args.outputs
+        )
     elif args.command == "port-forward":
         forward_ports()
     elif args.command == "ls":
diff --git a/rid/entrypoint/redim.py b/rid/entrypoint/redim.py
new file mode 100644
index 00000000..48d67dba
--- /dev/null
+++ b/rid/entrypoint/redim.py
@@ -0,0 +1,83 @@
+import json
+from pathlib import Path
+from typing import List, Union, Optional
+from rid.utils import load_json
+from copy import deepcopy
+import os
+
+from dflow import (
+    Workflow,
+    Step,
+    upload_artifact
+)
+
+from dflow.python import upload_packages
+from rid import SRC_ROOT
+upload_packages.append(SRC_ROOT)
+
+from rid.utils import normalize_resources
+from rid.superop.mcmc import MCMC
+from rid.op.mcmc_run import MCMCRun
+from rid.op.mcmc_plot import MCMCPlot
+
+def redim_rid(
+        rid_config: str,
+        machine_config: str,
+        models: Optional[Union[str, List[str]]] = None,
+        workflow_id_defined: Optional[str] = None
+    ):
+    with open(machine_config, "r") as mcg:
+        machine_config_dict = json.load(mcg)
+    resources = machine_config_dict["resources"]
+    tasks = machine_config_dict["tasks"]
+    normalized_resources = {}
+    for resource_type in resources.keys():
+        normalized_resources[resource_type] = normalize_resources(resources[resource_type])
+
+    mcmc_op = MCMC(
+        "mcmc",
+        MCMCRun,
+        MCMCPlot,
+        run_config = normalized_resources[tasks["mcmc_run_config"]],
+        plot_config = normalized_resources[tasks["mcmc_plot_config"]],
+        retry_times=1)
+    
+    jdata = deepcopy(load_json(rid_config))
+    mcmc_config = jdata["MCMC_Config"]
+    
+               
+    fe_models = []
+    assert isinstance(jdata["init_models"],list), "model input should be list."
+    for model in jdata["init_models"]:
+        fe_models.append(model)
+                
+    model_list = []
+            
+    if models is not None:
+        for model in models:
+            if os.path.basename(model) in fe_models:
+                model_list.append(model)
+
+    if len(model_list) == 0:
+        models_artifact = None
+    else:
+        models_artifact = upload_artifact([Path(p) for p in model_list], archive=None)
+    
+    task_names = []
+    for index in range(len(model_list)):
+        task_names.append("%03d"%index)
+
+    rid_steps = Step("rid-mcmc",
+            mcmc_op,
+            artifacts={
+                "models": models_artifact
+            },
+            parameters={
+                "mcmc_config": mcmc_config,
+                "task_names": task_names,
+                "block_tag" : "000"
+            },
+        )
+    wf = Workflow("rid-mcmc", pod_gc_strategy="OnPodSuccess", parallelism=10, id = workflow_id_defined)
+    wf.add(rid_steps)
+    wf.submit()
diff --git a/rid/entrypoint/reredim.py b/rid/entrypoint/reredim.py
new file mode 100644
index 00000000..b910bd0a
--- /dev/null
+++ b/rid/entrypoint/reredim.py
@@ -0,0 +1,109 @@
+import json
+from pathlib import Path
+from typing import List, Union, Optional
+from rid.utils import load_json
+from copy import deepcopy
+import os
+
+from dflow import (
+    Workflow,
+    Step,
+    upload_artifact
+)
+
+from dflow.python import upload_packages
+from rid import SRC_ROOT
+upload_packages.append(SRC_ROOT)
+
+from rid.utils import normalize_resources
+from rid.superop.mcmc import MCMC
+from rid.op.mcmc_run import MCMCRun
+from rid.op.mcmc_plot import MCMCPlot
+
+def reredim_rid(
+        workflow_id: str,
+        rid_config: str,
+        machine_config: str,
+        models: Optional[Union[str, List[str]]] = None,
+        workflow_id_defined: Optional[str] = None,
+        pod: Optional[str] = None
+    ):
+    with open(machine_config, "r") as mcg:
+        machine_config_dict = json.load(mcg)
+    resources = machine_config_dict["resources"]
+    tasks = machine_config_dict["tasks"]
+    normalized_resources = {}
+    for resource_type in resources.keys():
+        normalized_resources[resource_type] = normalize_resources(resources[resource_type])
+
+    mcmc_op = MCMC(
+        "mcmc",
+        MCMCRun,
+        MCMCPlot,
+        run_config = normalized_resources[tasks["mcmc_run_config"]],
+        plot_config = normalized_resources[tasks["mcmc_plot_config"]],
+        retry_times=None)
+    
+    jdata = deepcopy(load_json(rid_config))
+    mcmc_config = jdata["MCMC_Config"]
+    
+               
+    fe_models = []
+    assert isinstance(jdata["init_models"],list), "model input should be list."
+    for model in jdata["init_models"]:
+        fe_models.append(model)
+                
+    model_list = []
+            
+    if models is not None:
+        for model in models:
+            if os.path.basename(model) in fe_models:
+                model_list.append(model)
+
+    if len(model_list) == 0:
+        models_artifact = None
+    else:
+        models_artifact = upload_artifact([Path(p) for p in model_list], archive=None)
+    
+    task_names = []
+    for index in range(len(model_list)):
+        task_names.append("%03d"%index)
+
+    rid_steps = Step("rid-mcmc",
+            mcmc_op,
+            artifacts={
+                "models": models_artifact
+            },
+            parameters={
+                "mcmc_config": mcmc_config,
+                "task_names": task_names,
+                "block_tag" : "000"
+            },
+        )
+    
+    old_workflow = Workflow(id=workflow_id)
+    all_steps = old_workflow.query_step()
+
+    succeeded_steps = []
+    restart_flag = 1
+    for step in all_steps:
+        if step["type"] == "Pod":
+            pod_key = step["key"]
+            if pod_key is not None:
+                pod_key_list = pod_key.split("-")
+                pod_step = "-".join(pod_key_list[1:-1])
+                if pod is not None:
+                    if pod_step == pod:
+                        restart_flag = 0
+                else:
+                    if step["phase"] != "Succeeded":
+                        restart_flag = 0
+                    else:
+                        restart_flag = 1
+            
+            if restart_flag == 1:
+                succeeded_steps.append(step)
+                
+    wf = Workflow("rid-mcmc-continue", pod_gc_strategy="OnPodSuccess", parallelism=10, id = workflow_id_defined)
+    wf.add(rid_steps)
+    wf.submit(reuse_step=succeeded_steps)
diff --git a/rid/entrypoint/resubmit.py b/rid/entrypoint/resubmit.py
index f6fe7944..7438aad4 100644
--- a/rid/entrypoint/resubmit.py
+++ b/rid/entrypoint/resubmit.py
@@ -24,6 +24,7 @@ def resubmit_rid(
         topology: Optional[str],
         rid_config: str,
         machine_config: str,
+        workflow_id_defined: Optional[str] = None,
         iteration: Optional[str] = None,
         pod: Optional[str] = None,
         models: Optional[Union[str, List[str]]] = None,
@@ -50,6 +51,7 @@ def resubmit_rid(
         run_select_config = normalized_resources[tasks["run_select_config"]],
         prep_data_config = normalized_resources[tasks["prep_data_config"]],
         run_train_config = normalized_resources[tasks["run_train_config"]],
+        model_devi_config = normalized_resources[tasks["model_devi_config"]],
         workflow_steps_config = normalized_resources[tasks["workflow_steps_config"]],
         retry_times=None
     )
@@ -197,6 +199,6 @@ def resubmit_rid(
                 
                 if restart_flag == 1:
                     succeeded_steps.append(step)
-    wf = Workflow("reinforced-dynamics-continue", pod_gc_strategy="OnPodSuccess", parallelism=50)
+    wf = Workflow("reinforced-dynamics-continue", pod_gc_strategy="OnPodSuccess", parallelism=50, id = workflow_id_defined)
     wf.add(rid_steps)
     wf.submit(reuse_step=succeeded_steps)
diff --git a/rid/entrypoint/submit.py b/rid/entrypoint/submit.py
index 618993e1..2bdb98e3 100644
--- a/rid/entrypoint/submit.py
+++ b/rid/entrypoint/submit.py
@@ -29,6 +29,7 @@
 from rid.op.prep_data import CollectData, MergeData
 from rid.superop.blocks import IterBlock, InitBlock
 from rid.op.run_train import TrainModel
+from rid.op.run_model_devi import RunModelDevi
 from rid.op.adjust_trust_level import AdjustTrustLevel
 from rid.flow.loop import ReinforcedDynamics
 
@@ -42,6 +43,7 @@ def prep_rid_op(
     run_select_config,
     prep_data_config,
     run_train_config,
+    model_devi_config,
     workflow_steps_config,
     retry_times
     ):
@@ -85,7 +87,9 @@ def prep_rid_op(
         label_op,
         data_op,
         TrainModel,
+        RunModelDevi,
         run_train_config,
+        model_devi_config,
         retry_times=retry_times
     )
 
@@ -97,8 +101,10 @@ def prep_rid_op(
         data_op,
         AdjustTrustLevel,
         TrainModel,
+        RunModelDevi,
         workflow_steps_config,
         run_train_config,
+        model_devi_config,
         retry_times=retry_times)
     
     rid_op = ReinforcedDynamics(
@@ -114,6 +120,7 @@ def submit_rid(
         topology: Optional[str],
         rid_config: str,
         machine_config: str,
+        workflow_id_defined: Optional[str] = None,
         models: Optional[Union[str, List[str]]] = None,
         forcefield: Optional[str] = None,
         index_file: Optional[str] = None,
@@ -138,6 +145,7 @@ def submit_rid(
         run_select_config = normalized_resources[tasks["run_select_config"]],
         prep_data_config = normalized_resources[tasks["prep_data_config"]],
         run_train_config = normalized_resources[tasks["run_train_config"]],
+        model_devi_config = normalized_resources[tasks["model_devi_config"]],
         workflow_steps_config = normalized_resources[tasks["workflow_steps_config"]],
         retry_times=1
     )
@@ -255,6 +263,6 @@ def submit_rid(
             },
             parameters={}
         )
-    wf = Workflow("reinforced-dynamics", pod_gc_strategy="OnPodSuccess", parallelism=50)
+    wf = Workflow("reinforced-dynamics", pod_gc_strategy="OnPodSuccess", parallelism=50, id = workflow_id_defined)
     wf.add(rid_steps)
     wf.submit()
diff --git a/rid/flow/loop.py b/rid/flow/loop.py
index 52309abf..b7c6cabd 100644
--- a/rid/flow/loop.py
+++ b/rid/flow/loop.py
@@ -51,6 +51,7 @@ def __init__(
             "label_config": InputParameter(type=Dict),
             "tail": InputParameter(type=float, value=0.9),
             "train_config": InputParameter(type=Dict),
+            "type_map": InputParameter(type=List, value = []),
             "adjust_amplifier": InputParameter(type=float, value=1.5),
             "max_level_multiple": InputParameter(type=float, value=8.0),
         }
@@ -171,6 +172,7 @@ def _loop (
             "dt": steps.inputs.parameters["dt"],
             "output_freq": steps.inputs.parameters["output_freq"],
             "slice_mode": steps.inputs.parameters["slice_mode"],
+            "type_map": steps.inputs.parameters["type_map"],
             "label_config": steps.inputs.parameters["label_config"],
             "train_config": steps.inputs.parameters["train_config"]
         },
@@ -212,6 +214,7 @@ def _loop (
             "dt": steps.inputs.parameters["dt"],
             "output_freq": steps.inputs.parameters["output_freq"],
             "slice_mode": steps.inputs.parameters["slice_mode"],
+            "type_map": steps.inputs.parameters["type_map"],
             "label_config": steps.inputs.parameters["label_config"],
             "train_config": steps.inputs.parameters["train_config"]
         },
@@ -411,6 +414,7 @@ def _rid(
             "dt": prep_rid.outputs.parameters["dt"],
             "output_freq": prep_rid.outputs.parameters["output_freq"],
             "slice_mode": prep_rid.outputs.parameters["slice_mode"],
+            "type_map": prep_rid.outputs.parameters["type_map"],
             "label_config": prep_rid.outputs.parameters["label_config"],
             "train_config": prep_rid.outputs.parameters["train_config"]
         },
@@ -457,6 +461,7 @@ def _rid(
             "dt": prep_rid.outputs.parameters["dt"],
             "output_freq": prep_rid.outputs.parameters["output_freq"],
             "slice_mode": prep_rid.outputs.parameters["slice_mode"],
+            "type_map": prep_rid.outputs.parameters["type_map"],
             "label_config": prep_rid.outputs.parameters["label_config"],
             "train_config": prep_rid.outputs.parameters["train_config"]
         },
diff --git a/rid/mcmc/__init__.py b/rid/mcmc/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/rid/mcmc/walker.py b/rid/mcmc/walker.py
new file mode 100644
index 00000000..4684706d
--- /dev/null
+++ b/rid/mcmc/walker.py
@@ -0,0 +1,127 @@
+import numpy as np
+import tensorflow as tf
+
+# kinetic enery in eV
+kbT = (8.617343E-5) * 300 
+beta = 1.0 / kbT
+
+# ev to kj/mol
+f_cvt = 96.485
+
+def distance(cv1,cv2):
+    d = np.sum((cv1 - cv2)**2,axis=1)
+    return d
+
+def s(cv, lm, xlist):
+    s = 0
+    n = 1e-10
+    for i in range(1,10):
+       s +=  i*np.exp(-lm*distance(cv,xlist[i-1]))
+       n += np.exp(-lm*distance(cv,xlist[i-1]))
+    return s/n
+
+def z(cv, lm, xlist):
+    z = 0
+    for i in range(9):
+        z += np.exp(-lm*distance(cv,xlist[i]))
+    p = np.ones(shape=z.shape)*1e-10
+    z = np.log(np.max((z,p),axis=0))
+    return -z/lm
+
+class Walker(object):
+    def __init__(self, fd, nw, sess, type, cv_lower, cv_upper):
+        self.type = type
+        self._full_dim = fd
+        self._num_walker = nw
+        self._move_scale = 0.5
+        self._sess = sess
+        self._shape = (self._num_walker, self._full_dim)
+        self._cv_range = cv_upper
+        self._cv_shift = cv_lower
+        # absolute coordinate
+        self._position = np.random.uniform(size=self._shape)*(np.array(self._cv_range) - np.array(self._cv_shift)) + np.array(self._cv_shift)
+        self._energy = np.zeros([self._num_walker])
+        self._force = np.zeros([self._num_walker, self._full_dim])
+        self._sample_step = 20
+        self._acp_ratio_lb = 0.15
+        self._acp_ratio_ub = 0.75
+        self.max_scale = 5
+        self.min_scale = 0.01
+        self.inc_scale_fac = 1.25
+
+    def sample(self, compute_ef, inter_step=1):
+        acp_ratio = []
+        self._energy, self._force = compute_ef(self._sess, self._position)
+        for _ in range(inter_step):
+            if self.type == "dih":
+                position_new = np.mod(self._position + np.random.normal(scale=self._move_scale,size=self._shape), 2*np.pi)
+            elif self.type == "dis":
+                position_new = self._position + np.random.normal(scale=self._move_scale*0.08, size=self._shape)
+                indices = np.where(position_new < 0.1)
+                position_new[indices]=0.13 + np.random.normal(scale=0.01, size=indices[0].shape)
+                indices = np.where(position_new > 9.9)
+                position_new[indices]=9.87 + np.random.normal(scale=0.01, size=indices[0].shape)
+            else:
+                raise ValueError("Undefined cv type, only support 'dih' and 'dis' type")
+            energy_new, force_new = compute_ef(self._sess, position_new)
+            # in case of overflow
+            prob_ratio = np.exp(np.minimum(- beta * (energy_new - self._energy), 0))
+            idx = np.random.uniform(size=self._num_walker) < np.reshape(prob_ratio, [-1])
+
+            self._position[idx, :] = position_new[idx, :]
+            self._energy[idx] = energy_new[idx]
+            self._force[idx] = force_new[idx]
+            acp_ratio.append(np.mean(idx))
+            # elif self.type == "dis":
+            #     for ii in range(self._num_walker):
+            #         if (position_new[ii,:] > self._cv_shift).all() and (position_new[ii,:] < self._cv_range).all() and idx[ii]:
+            #             self._position[ii, :] = position_new[ii, :]
+            #             self._energy[ii] = energy_new[ii]
+            #             self._force[ii] = force_new[ii]
+            #             acp_ratio.append(np.mean(idx))
+            # else:
+            #     raise ValueError("Undefined cv type, only support 'dih' and 'dis' type")
+        acp_ratio = np.mean(acp_ratio)
+        if acp_ratio > self._acp_ratio_ub:
+            # move_scale is too small
+            self._move_scale = np.min((self._move_scale*self.inc_scale_fac,self.max_scale),axis=0)
+            print(
+                "Increase move_scale to %s due to high acceptance ratio: %f" % (
+                    self._move_scale, acp_ratio))
+            # print(self._position[:5, :, :])
+        elif acp_ratio < self._acp_ratio_lb:
+            # move_scale is too large
+            self._move_scale = np.max((self._move_scale/self.inc_scale_fac,self.min_scale),axis=0)
+            print(
+                "Decrease move_scale to %s due to low acceptance ratio: %f" % (
+                    self._move_scale, acp_ratio))
+        return self._position, self._energy, self._force
+
+# project on all 1d points
+def my_hist1d(pp, xx, delta, fd):
+    my_hist = np.zeros((fd, len(xx)))
+    for ii in range(pp.shape[0]):   ###trj_num
+        for jj in range(fd):        ###cv_num
+            my_hist[jj, np.int(pp[ii,jj]//delta)] += 1
+    my_hist /= (pp.shape[0] * delta)
+    return my_hist
+
+# project on specific 2d point
+def my_hist2d(pp, xx, yy, delta, cv1, cv2):
+    my_hist = np.zeros((1, len(xx), len(yy)))
+    for ii in range(pp.shape[0]):
+        my_hist[0, np.int(pp[ii,cv1]//delta), np.int(pp[ii,cv2]//delta)] += 1
+    my_hist /= (pp.shape[0] * delta * delta)
+    return my_hist
+
+# project on path CV
+def my_hist2d_path(pp, xx, yy, delta, lm, xlist, proj_index):
+    my_hist = np.zeros((1, len(xx), len(yy)))
+    cv = pp[:,proj_index]
+    slist = s(cv, lm, xlist)
+    zlist = z(cv, lm, xlist)
+    for ii in range(cv.shape[0]):
+        if 0 < slist[ii] < 10 and -10 < zlist[ii] < 0:
+            my_hist[0, np.int((slist[ii])//delta), np.int((zlist[ii]+10)//delta)] += 1
+    my_hist /= (pp.shape[0] * delta * delta)
+    return my_hist
\ No newline at end of file
diff --git a/rid/op/mcmc_plot.py b/rid/op/mcmc_plot.py
new file mode 100644
index 00000000..ce7a55d9
--- /dev/null
+++ b/rid/op/mcmc_plot.py
@@ -0,0 +1,155 @@
+from typing import List, Optional, Dict
+from pathlib import Path
+import numpy as np
+from dflow.python import (
+    OP,
+    OPIO,
+    OPIOSign,
+    Artifact,
+    Parameter,
+    BigParameter
+)
+from rid.utils import save_txt, set_directory
+from rid.constants import mcmc_1cv_dir_name, mcmc_1cv_name, mcmc_2cv_name,mcmc_1cv_fig, mcmc_2cv_fig, mcmc_2cv_fig_separate
+from matplotlib import pyplot as plt
+import os
+
+class MCMCPlot(OP):
+    """
+    `MCMC_Plot` plot the reduced free energy surface produced by MCMC run.
+    """
+
+    @classmethod
+    def get_input_sign(cls):
+        return OPIOSign(
+            {
+                "mcmc_1cv": Artifact(List[Path]),
+                "mcmc_2cv": Artifact(List[Path]),
+                "mcmc_config": BigParameter(Dict)
+            }
+        )
+
+    @classmethod
+    def get_output_sign(cls):
+        return OPIOSign(
+            {
+                "mcmc_fig": Artifact(Path, archive = None)
+            }
+        )
+
+    @OP.exec_sign_check
+    def execute(
+        self,
+        op_in: OPIO,
+    ) -> OPIO:
+
+        r"""Execute the OP.
+        
+        Parameters
+        ----------
+        op_in : dict
+            Input dict with components:
+
+            - `mcmc_config`: Dict,
+            - `mcmc_cv`: Artifact(List[Path]),
+          
+        Returns
+        -------
+            Output dict with components:
+        
+            - `mcmc_fig`: (`Artifact(Path)`)
+        """
+        mcmc_config = op_in["mcmc_config"]
+        cv_type = mcmc_config["cv_type"]
+        bins = mcmc_config["bins"]
+        cv_dim = mcmc_config["cv_dimension"]
+        
+        task_path = Path("mcmc_fig")
+        task_path.mkdir(exist_ok=True, parents=True)
+        with set_directory(task_path):
+            fourdata_1cv = []
+            for dir in op_in["mcmc_1cv"]:
+                all_1cv = []
+                for cv_index in range(cv_dim):
+                    file = dir/mcmc_1cv_name.format(tag = cv_index)
+                    pmf = np.loadtxt(file)
+                    pmf = pmf - np.min(pmf)
+                    proj_iterations = pmf.shape[0]
+                    all_1cv.append(pmf)
+                fourdata_1cv.append(all_1cv)
+            avedata_1cv = np.mean(np.array(fourdata_1cv),axis=0)
+            print("avedata shape", avedata_1cv.shape)
+            # make 1cv plot
+            if not os.path.exists(mcmc_1cv_dir_name):
+                os.makedirs(mcmc_1cv_dir_name)
+            if cv_type == "dih":
+                xedges = np.linspace(0, 2*np.pi, bins)
+            elif cv_type == "dis":
+                xedges = np.linspace(0, 10, bins)
+            for cv_index in range(cv_dim):
+                plt.figure(figsize=(8, 6))
+                for proj_iter in range(proj_iterations):
+                    print("xedges shape", xedges.shape)
+                    plt.plot(xedges,avedata_1cv[cv_index][proj_iter], label = proj_iter)
+                plt.xlabel(r'cv')
+                plt.ylabel(r'free energy (kcal/mol)')
+                plt.legend()
+                plt.savefig(mcmc_1cv_dir_name+"/"+mcmc_1cv_fig.format(tag = cv_index),dpi=600,bbox_inches='tight')
+                
+            fourdata_2cv = []
+            for file in op_in["mcmc_2cv"]:
+                pmf = np.loadtxt(file)
+                pmf = pmf - np.min(pmf)
+                allda=[]
+                if cv_type == "dih":
+                    xedges = np.linspace(0, 2*np.pi, bins)
+                    yedges = np.linspace(0, 2*np.pi, bins)
+                elif cv_type == "dis":
+                    xedges = np.linspace(0, 10, bins)
+                    yedges = np.linspace(0, 10, bins)
+                for i in range(bins):
+                    for j in range(bins):
+                        temp = []
+                        temp.append(xedges[i])
+                        temp.append(yedges[j])
+                        temp.append(pmf[i][j])
+                        temp = np.array(temp)
+                        allda.append(temp)
+                newarray=np.array(allda)
+                idex=np.lexsort([newarray[:,1], newarray[:,0]])
+                sorted_data = newarray[idex,:]
+                fourdata_2cv.append(sorted_data[:,2])
+            avedata_2cv = np.mean(np.array(fourdata_2cv),axis=0)    
+            # make 2cv plot
+            if cv_type == "dih":
+                xedges = np.linspace(0, 2*np.pi, bins)
+                yedges = np.linspace(0, 2*np.pi, bins)
+            elif cv_type == "dis":
+                xedges = np.linspace(0, 10, bins)
+                yedges = np.linspace(0, 10, bins)
+            plt.figure(figsize=(8, 6))
+            cmap = plt.cm.get_cmap("jet_r")
+            CS = plt.contourf(xedges,yedges,avedata_2cv.reshape(bins,bins),levels = np.linspace(0,6,61),cmap=cmap,extend="max")
+            cbar = plt.colorbar(CS)
+            cbar.ax.tick_params(labelsize=8) 
+            cbar.ax.set_title('kcal/mol',fontsize=8)
+            plt.xlabel(r'CV index 1')
+            plt.ylabel(r'CV index 2')
+            plt.savefig(mcmc_2cv_fig,dpi=600,bbox_inches='tight')
+            for iii in range(len(fourdata_2cv)):
+                fig = plt.figure(figsize=(8, 6))
+                cmap = plt.cm.get_cmap("jet_r")
+                CS = plt.contourf(xedges,yedges,fourdata_2cv[iii].reshape(bins,bins),levels = np.linspace(0,6,61),cmap=cmap,extend="max")
+                cbar = plt.colorbar(CS)
+                cbar.ax.tick_params(labelsize=8) 
+                cbar.ax.set_title('kcal/mol',fontsize=8)
+                plt.xlabel(r'CV index 1')
+                plt.ylabel(r'CV index 2')
+                plt.savefig(mcmc_2cv_fig_separate.format(tag=iii),dpi=600,bbox_inches='tight')
+                
+        op_out = OPIO(
+            {
+               "mcmc_fig": task_path
+            }
+        )
+        return op_out
\ No newline at end of file
diff --git a/rid/op/mcmc_run.py b/rid/op/mcmc_run.py
new file mode 100644
index 00000000..198ef6e7
--- /dev/null
+++ b/rid/op/mcmc_run.py
@@ -0,0 +1,172 @@
+from typing import List, Optional, Dict
+from pathlib import Path
+import numpy as np
+from dflow.python import (
+    OP,
+    OPIO,
+    OPIOSign,
+    Artifact,
+    Parameter,
+    BigParameter
+)
+from rid.utils import save_txt, set_directory
+from rid.mcmc.walker import Walker, my_hist1d, my_hist2d, my_hist2d_path
+from rid.select.model_devi import test_ef
+from rid.common.tensorflow.graph import load_graph
+try:
+    import tensorflow.compat.v1 as tf
+    tf.disable_v2_behavior()
+except ImportError:
+    import tensorflow as tf
+from rid.constants import mcmc_1cv_name, mcmc_1cv_dir_name, mcmc_2cv_name, kb, f_cvt, kcal2kj
+import os
+
+class MCMCRun(OP):
+    """
+    `MCMC_Run` performs Markov chain Monte carlo sampling to reduce the dimension of Free energy surface (represented by neural network)
+    """
+
+    @classmethod
+    def get_input_sign(cls):
+        return OPIOSign(
+            {
+                "task_names": str,
+                "mcmc_config": BigParameter(Dict),
+                "models": Artifact(Path, archive = None)
+            }
+        )
+
+    @classmethod
+    def get_output_sign(cls):
+        return OPIOSign(
+            {
+                "mcmc_1cv": Artifact(Path, archive = None),
+                "mcmc_2cv": Artifact(Path, archive = None)
+            }
+        )
+
+    @OP.exec_sign_check
+    def execute(
+        self,
+        op_in: OPIO,
+    ) -> OPIO:
+
+        r"""Execute the OP.
+        
+        Parameters
+        ----------
+        op_in : dict
+            Input dict with components:
+
+            - `task_name`: str,
+            - `models`: Artifact(List[Path], optional=True),
+          
+        Returns
+        -------
+            Output dict with components:
+        
+            - `fes_fig`: (`Artifact(Path)`)
+        """
+        graph = load_graph(op_in["models"])
+        mcmc_config = op_in["mcmc_config"]
+        temperature = mcmc_config["temperature"]
+        fd = mcmc_config["cv_dimension"]
+        ns = mcmc_config["numb_steps"]
+        nw = mcmc_config["numb_walkers"]
+        # kinetic enery in eV
+        kbT = kb * temperature
+        beta = 1.0 / kbT
+        if "cv_upper_bound" in mcmc_config:
+            cv_upper = mcmc_config["cv_upper_bound"]
+        else:
+            cv_upper = 2*np.pi
+        if "cv_lower_bound" in mcmc_config:
+            cv_lower = mcmc_config["cv_lower_bound"]
+        else:
+            cv_lower = 0
+        proj_info = mcmc_config["proj_info"]
+        proj_mode = proj_info["proj_mode"]
+        proj_cv_index = proj_info["proj_cv_index"]
+        if "path_list" in proj_info:
+            path_list = np.array(proj_info["path_list"])
+            path_lm = proj_info["path_lm"]
+        cv_type = mcmc_config["cv_type"]
+        bins = mcmc_config["bins"]
+        
+        if cv_type == "dih":
+            xx = np.linspace(0,2* np.pi, bins)
+            yy = np.linspace(0,2* np.pi, bins)
+            pp_hist = np.zeros((fd, len(xx)))
+            pp_hist2d = np.zeros((1, len(xx), len(yy)))
+            delta = 2.0 * np.pi / (bins-1)
+        elif cv_type == "dis":
+            xx = np.linspace(0,10, bins)
+            yy = np.linspace(0,10, bins)
+            pp_hist = np.zeros((fd, len(xx)))
+            pp_hist2d = np.zeros((1, len(xx), len(yy)))
+            delta = 10.0 / (bins-1)
+        else:
+            raise ValueError("Undefined cv type, only support 'dih' and 'dis' type")
+        
+        task_path = Path(op_in["task_names"])
+        task_path.mkdir(exist_ok=True, parents=True)
+
+        mcmc_2cv_path = None
+        with set_directory(task_path):
+            with tf.Session(graph = graph) as sess:        
+                walker = Walker(fd, nw, sess, cv_type, cv_lower=cv_lower, cv_upper=cv_upper)
+                for _ in range(10000):
+                    pp, ee, ff = walker.sample(test_ef)
+
+                for ii in range(ns+1):
+                    pp, ee, ff = walker.sample(test_ef)
+                    pp_hist_new = my_hist1d(pp, xx, delta, fd)
+                    pp_hist = (pp_hist * ii + pp_hist_new) / (ii+1)
+                    if not os.path.exists(mcmc_1cv_dir_name):
+                        os.makedirs(mcmc_1cv_dir_name)
+                    if np.mod(ii,int(ns/5)) == 0:
+                        zz = -np.log(pp_hist+1e-7)/beta
+                        # convert ev to kcal/mol
+                        zz *= f_cvt/kcal2kj
+                        zz = zz - np.min(zz)      
+                        for jj in range(fd):
+                            fp = open(mcmc_1cv_dir_name+"/"+mcmc_1cv_name.format(tag=jj), "a")
+                            for temp in zz[jj]:
+                                fp.write(str(temp)+'    ')
+                            fp.write('\n')
+                            fp.close()
+                        
+                    if proj_mode == "cv":
+                        if len(proj_cv_index) == 2:
+                            cv1 = proj_cv_index[0]
+                            cv2 = proj_cv_index[1]
+                            ##certain 2d
+                            pp_hist_new2d = my_hist2d(pp, xx, yy, delta, cv1, cv2)
+                            pp_hist2d = (pp_hist2d * ii + pp_hist_new2d) / (ii+1)
+                            if ii == ns:
+                                zz2d = np.transpose(-np.log(pp_hist2d+1e-10), (0,2,1))/beta
+                                # convert ev to kcal/mol
+                                zz2d *= f_cvt/kcal2kj
+                                zz2d = zz2d - np.min(zz2d)
+                                np.savetxt(mcmc_2cv_name,zz2d[0])
+                    elif proj_mode == "path":
+                        pp_hist_new2d_path = my_hist2d_path(pp, xx, yy, delta, path_lm, path_list, proj_cv_index)
+                        pp_hist2d = (pp_hist2d * ii + pp_hist_new2d_path) / (ii+1)
+                        if ii == ns:
+                            zz2d = np.transpose(-np.log(pp_hist2d+1e-10), (0,2,1))/beta
+                            # convert ev to kcal/mol
+                            zz2d *= f_cvt/kcal2kj
+                            zz2d = zz2d - np.min(zz2d)
+                            np.savetxt(mcmc_2cv_name,zz2d[0])
+                            
+        if os.path.exists(task_path.joinpath(mcmc_2cv_name)):
+            mcmc_2cv_path = task_path.joinpath(mcmc_2cv_name)
+            
+            
+        op_out = OPIO(
+            {
+               "mcmc_1cv": task_path.joinpath(mcmc_1cv_dir_name),
+               "mcmc_2cv": mcmc_2cv_path
+            }
+        )
+        return op_out
\ No newline at end of file
diff --git a/rid/op/prep_data.py b/rid/op/prep_data.py
index 80c7d7f0..9b314764 100644
--- a/rid/op/prep_data.py
+++ b/rid/op/prep_data.py
@@ -43,10 +43,11 @@ def execute(
         ) -> OPIO:
         cv_forces = []
         for idx in range(len(op_in["cv_forces"])):
-            cv_force = load_txt(op_in["cv_forces"][idx])
-            cv_forces = np.append(cv_forces, cv_force)
+            if op_in["cv_forces"][idx]:
+                cv_force = load_txt(op_in["cv_forces"][idx])
+                cv_forces = np.append(cv_forces, cv_force)
         if op_in["cv_forces"]:
-            cv_forces = np.reshape(cv_forces, [len(op_in["cv_forces"]), -1])
+            cv_forces = np.reshape(cv_forces, [-1, len(cv_force)])
             data = cv_forces
         else:
             data = np.array([])
diff --git a/rid/op/prep_exploration.py b/rid/op/prep_exploration.py
index a749614c..f2f9c0f9 100644
--- a/rid/op/prep_exploration.py
+++ b/rid/op/prep_exploration.py
@@ -3,7 +3,8 @@
     OP,
     OPIO,
     OPIOSign,
-    Artifact
+    Artifact,
+    BigParameter
 )
 
 from typing import List, Dict, Union
@@ -51,8 +52,8 @@ def get_input_sign(cls):
                 "cv_file": Artifact(List[Path], optional=True),
                 "trust_lvl_1": float,
                 "trust_lvl_2": float,
-                "exploration_config": Dict,
-                "cv_config": Dict,
+                "exploration_config": BigParameter(Dict),
+                "cv_config": BigParameter(Dict),
                 "task_name": str,
                 "block_tag": str
             }
@@ -62,7 +63,7 @@ def get_input_sign(cls):
     def get_output_sign(cls):
         return OPIOSign(
             {
-                "task_path": Artifact(Path),
+                "task_path": Artifact(Path, archive = None),
                 "cv_dim": int
             }
         )
diff --git a/rid/op/prep_label.py b/rid/op/prep_label.py
index 33f52f41..dfc99b9c 100644
--- a/rid/op/prep_label.py
+++ b/rid/op/prep_label.py
@@ -3,7 +3,8 @@
     OPIO,
     OPIOSign,
     Artifact,
-    Parameter
+    Parameter,
+    BigParameter
 )
 
 import json
@@ -29,7 +30,7 @@ def get_input_sign(cls):
         return OPIOSign(
             {
                 "confs": Artifact(List[Path], optional=True),
-                "conf_tags": Parameter(List, default=[])
+                "conf_tags": Artifact(List[Path], optional=True)
             }
         )
 
@@ -38,7 +39,7 @@ def get_output_sign(cls):
         return OPIOSign(
             {
                 "if_continue": int,
-                "conf_tags": List
+                "conf_tags": BigParameter(List)
             }
         )
 
@@ -69,10 +70,9 @@ def execute(
 
             tags = {}
             for tag in op_in["conf_tags"]:
-                if isinstance(tag, Dict):
-                    tags.update(tag)
-                elif isinstance(tag, str):
-                    tags.update(json.loads(tag))
+                if isinstance(tag,Path):
+                    with open(tag,"r") as f:
+                        tags.update(json.load(f))
                 else:
                     raise RuntimeError("Unkown Error.")
             
@@ -104,9 +104,9 @@ def get_input_sign(cls):
                 "topology": Artifact(Path, optional=True),
                 "conf": Artifact(Path),
                 "cv_file": Artifact(List[Path], optional=True),
-                "label_config": Dict,
-                "cv_config": Dict,
-                "task_name": str,
+                "label_config": BigParameter(Dict),
+                "cv_config": BigParameter(Dict),
+                "task_name": BigParameter(str),
                 "at": Artifact(Path, optional=True)
             }
         )
@@ -115,7 +115,7 @@ def get_input_sign(cls):
     def get_output_sign(cls):
         return OPIOSign(
             {
-                "task_path": Artifact(Path),
+                "task_path": Artifact(Path, archive = None),
             }
         )
 
diff --git a/rid/op/prep_rid.py b/rid/op/prep_rid.py
index 320d7b2a..de9a7319 100644
--- a/rid/op/prep_rid.py
+++ b/rid/op/prep_rid.py
@@ -6,10 +6,11 @@
     OP,
     OPIO,
     OPIOSign,
-    Artifact
+    Artifact,
+    BigParameter
 )
 from rid.utils import load_json
-from rid.constants import model_tag_fmt, init_conf_name, init_input_name, walker_tag_fmt
+from rid.constants import model_tag_fmt, init_conf_gmx_name, init_conf_lmp_name,init_input_name, walker_tag_fmt
 
 
 logging.basicConfig(
@@ -21,22 +22,39 @@
 logger = logging.getLogger(__name__)
 
 
-def prep_confs(confs, numb_walkers):
+def prep_confs(confs, numb_walkers,sampler_type):
     numb_confs = len(confs)
     conf_list = []
-    if numb_confs < numb_walkers:
-        logger.info("Number of confs is smaller than number of walkers. Copy replicas up to number of walkers.")
+    if sampler_type == "gmx":
+        if numb_confs < numb_walkers:
+            logger.info("Number of confs is smaller than number of walkers. Copy replicas up to number of walkers.")
+            for idx in range(numb_walkers):
+                shutil.copyfile(confs[idx%numb_confs], init_conf_gmx_name.format(idx=idx))
+        elif numb_confs > numb_walkers:
+            logger.info("Number of confs is greater than number of walkers. Only use the fist `numb_walkers` confs.")
+            for idx in range(numb_walkers):
+                shutil.copyfile(confs[idx], init_conf_gmx_name.format(idx=idx))
+        else:
+            for idx in range(numb_walkers):
+                shutil.copyfile(confs[idx], init_conf_gmx_name.format(idx=idx))
         for idx in range(numb_walkers):
-            shutil.copyfile(confs[idx%numb_confs], init_conf_name.format(idx=idx))
-    elif numb_confs > numb_walkers:
-        logger.info("Number of confs is greater than number of walkers. Only use the fist `numb_walkers` confs.")
+            conf_list.append(Path(init_conf_gmx_name.format(idx=idx)))
+    elif sampler_type == "lmp":
+        if numb_confs < numb_walkers:
+            logger.info("Number of confs is smaller than number of walkers. Copy replicas up to number of walkers.")
+            for idx in range(numb_walkers):
+                shutil.copyfile(confs[idx%numb_confs], init_conf_lmp_name.format(idx=idx))
+        elif numb_confs > numb_walkers:
+            logger.info("Number of confs is greater than number of walkers. Only use the fist `numb_walkers` confs.")
+            for idx in range(numb_walkers):
+                shutil.copyfile(confs[idx], init_conf_lmp_name.format(idx=idx))
+        else:
+            for idx in range(numb_walkers):
+                shutil.copyfile(confs[idx], init_conf_lmp_name.format(idx=idx))
         for idx in range(numb_walkers):
-            shutil.copyfile(confs[idx], init_conf_name.format(idx=idx))
+            conf_list.append(Path(init_conf_lmp_name.format(idx=idx)))
     else:
-        for idx in range(numb_walkers):
-            shutil.copyfile(confs[idx], init_conf_name.format(idx=idx))
-    for idx in range(numb_walkers):
-        conf_list.append(Path(init_conf_name.format(idx=idx)))
+        raise ValueError("Invalid sampler type, only support gmx and lmp")
     return conf_list
 
 
@@ -65,12 +83,12 @@ def get_output_sign(cls):
                 "numb_iters": int,
                 "numb_walkers": int,
                 "numb_models": int,
-                "confs": Artifact(List[Path]),
+                "confs": Artifact(List[Path],archive = None),
                 "walker_tags": List,
                 "model_tags": List,
                 
-                "exploration_config": Dict,
-                "cv_config": Dict,
+                "exploration_config": BigParameter(Dict),
+                "cv_config": BigParameter(Dict),
                 "trust_lvl_1": List[float],
                 "trust_lvl_2": List[float],
                 "cluster_threshold": List[float],
@@ -83,8 +101,9 @@ def get_output_sign(cls):
                 "dt": float,
                 "output_freq": float,
                 "slice_mode": str,
-                "label_config": Dict,
-                "train_config": Dict
+                "type_map": List,
+                "label_config": BigParameter(Dict),
+                "train_config": BigParameter(Dict)
             }
         )
 
@@ -142,7 +161,10 @@ def execute(
         train_config = jdata.pop("Train")
         numb_models = train_config.pop("numb_models")
         numb_iters = jdata.pop("numb_iters")
-        conf_list = prep_confs(op_in["confs"], numb_walkers)
+        exploration_config = jdata.pop("ExploreMDConfig")
+        
+        sampler_type = exploration_config["type"]
+        conf_list = prep_confs(op_in["confs"], numb_walkers, sampler_type)
 
         walker_tags = []
         model_tags = []
@@ -151,7 +173,6 @@ def execute(
         for idx in range(numb_models):
             model_tags.append(model_tag_fmt.format(idx=idx))
         
-        exploration_config = jdata.pop("ExploreMDConfig")
         dt = exploration_config["dt"]
         output_freq = exploration_config["output_freq"]
         cv_config = jdata.pop("CV")
@@ -170,6 +191,11 @@ def execute(
         cluster_threshold = selection_config.pop("cluster_threshold")
         cluster_threshold_list = [cluster_threshold for _ in range(numb_walkers)]
         
+        if "type_map" in selection_config:
+            type_map = selection_config["type_map"]
+        else:
+            type_map = []
+        
         op_out = OPIO(
             {
                 "numb_iters": numb_iters,
@@ -193,6 +219,7 @@ def execute(
                 "dt": dt,
                 "output_freq": output_freq,
                 "slice_mode": selection_config.pop("slice_mode"),
+                "type_map": type_map,
                 "label_config": label_config,
                 "train_config": train_config
             }
diff --git a/rid/op/prep_select.py b/rid/op/prep_select.py
index 88c058bd..8a9f4ed1 100644
--- a/rid/op/prep_select.py
+++ b/rid/op/prep_select.py
@@ -53,9 +53,9 @@ def get_output_sign(cls):
             {
                 "numb_cluster": int,
                 "cluster_threshold": float,
-                "cluster_fig": Artifact(Path),
-                "cluster_selection_index": Artifact(Path),
-                "cluster_selection_data": Artifact(Path)
+                "cluster_fig": Artifact(Path, archive = None),
+                "cluster_selection_index": Artifact(Path, archive = None),
+                "cluster_selection_data": Artifact(Path, archive = None)
             }
         )
 
diff --git a/rid/op/run_exploration.py b/rid/op/run_exploration.py
index 65df3a07..0193b660 100644
--- a/rid/op/run_exploration.py
+++ b/rid/op/run_exploration.py
@@ -8,7 +8,8 @@
     OP,
     OPIO,
     OPIOSign,
-    Artifact
+    Artifact,
+    BigParameter
 )
 from rid.constants import (
         gmx_conf_name,
@@ -24,12 +25,17 @@
         lmp_conf_out,
         bias_fig,
         model_devi_fig,
+        dp_model_devi_fig,
+        dp_sel_ndx,
+        sel_lmp_name,
         lmp_input_name
     )
 from rid.utils import run_command, set_directory, list_to_string
 from rid.common.lammps.command import final_dump
 from rid.common.sampler.command import get_grompp_cmd, get_mdrun_cmd
+from rid.common.mol_dpdata import slice_dump
 from rid.select.model_devi import make_std
+from rid.mcmc.walker import s,z
 
 
 logging.basicConfig(
@@ -55,7 +61,7 @@ def get_input_sign(cls):
             {
                 "task_path": Artifact(Path),
                 "forcefield": Artifact(Path, optional=True),
-                "exploration_config": Dict,
+                "exploration_config": BigParameter(Dict),
                 "models": Artifact(List[Path], optional=True),
                 "index_file": Artifact(Path, optional=True),
                 "dp_files": Artifact(List[Path], optional=True),
@@ -68,12 +74,17 @@ def get_input_sign(cls):
     def get_output_sign(cls):
         return OPIOSign(
             {
-                "plm_out": Artifact(Path),
-                "bias_fig": Artifact(Path, optional=True),
-                "model_devi_fig": Artifact(Path,optional=True),
-                "md_log": Artifact(Path),
-                "trajectory": Artifact(Path),
-                "conf_out": Artifact(Path)
+                "plm_out": Artifact(Path, archive = None),
+                "bias_fig": Artifact(Path, optional=True, archive = None),
+                "model_devi_fig": Artifact(Path,optional=True, archive = None),
+                "dp_model_devi_fig": Artifact(Path,optional=True, archive = None),
+                "dp_model_devi": Artifact(Path,optional=True, archive = None),
+                "dp_selected_indices": Artifact(Path,optional=True, archive = None),
+                "dp_selected_confs": Artifact(List[Path], archive = None),
+                "projected_fig": Artifact(Path,optional=True, archive = None),
+                "md_log": Artifact(Path, archive = None),
+                "trajectory": Artifact(Path, archive = None),
+                "conf_out": Artifact(Path, archive = None)
             }
         )
 
@@ -159,24 +170,28 @@ def execute(
                 for file in op_in["dp_files"]:
                     if not os.path.islink(file.name):
                         os.symlink(file, file.name)
-            if op_in["index_file"] is not None:
-                if not os.path.islink(op_in["index_file"].name):
-                    os.symlink(op_in["index_file"], op_in["index_file"].name)
+
             if op_in["cv_file"] is not None:
                 for file in op_in["cv_file"]:
                     if file.name != "colvar":
                         if not os.path.islink(file.name):
                             os.symlink(file, file.name)
-            if op_in["inputfile"] is not None:
-                for file in op_in["inputfile"]:
-                    if file.name == inputfile_name:
-                        if not os.path.islink(file.name):
-                            os.symlink(file, file.name)
+                            
+            if op_in["exploration_config"]["type"] == "lmp":
+                if op_in["inputfile"] is not None:
+                    for file in op_in["inputfile"]:
+                        if file.name == inputfile_name:
+                            if not os.path.islink(file.name):
+                                os.symlink(file, file.name)
             
-            if op_in["dp_files"] is not None:
-                for dp_file in op_in["dp_files"]:
-                    if (dp_file.name).endswith("json"):
-                        os.environ["GMX_DEEPMD_INPUT_JSON"] = dp_file.name
+            if op_in["exploration_config"]["type"] == "gmx":
+                if op_in["index_file"] is not None:
+                    if not os.path.islink(op_in["index_file"].name):
+                        os.symlink(op_in["index_file"], op_in["index_file"].name)
+                if op_in["dp_files"] is not None:
+                    for dp_file in op_in["dp_files"]:
+                        if (dp_file.name).endswith("json"):
+                            os.environ["GMX_DEEPMD_INPUT_JSON"] = dp_file.name
                         
             if grompp_cmd is not None:
                 logger.info(list_to_string(grompp_cmd, " "))
@@ -193,6 +208,77 @@ def execute(
             
             bias_fig_file = None
             model_devi_fig_file = None
+            
+            dp_model_devi_fig_file = None
+            dp_model_devi_file = None
+            dp_sel_ndx_file = None
+            dp_conf_list = []
+            if "dp_model_devi_out" in op_in["exploration_config"]:
+                # plot DP potential model deviation
+                dp_model_devi_file_name = op_in["exploration_config"]["dp_model_devi_out"]
+                data = np.loadtxt(dp_model_devi_file_name)
+                f_list = data[:,4]
+                cv_data = np.loadtxt("plm.out")[:,2:]
+                assert (len(cv_data) == len(f_list))
+                trust_lvl_hi = op_in["exploration_config"]["dp_model_devi_f_hi"]
+                trust_lvl_lo = op_in["exploration_config"]["dp_model_devi_f_lo"]
+                max_selection = op_in["exploration_config"]["dp_model_max_selection"]
+                high_devi = np.where((f_list>trust_lvl_lo) & (f_list < trust_lvl_hi))
+                
+                if len(high_devi[0]) > max_selection:
+                    high_devi = np.random.choice(high_devi[0],max_selection,replace=False)
+                else:
+                    high_devi = high_devi[0]
+                np.savetxt(dp_sel_ndx,high_devi,fmt="%d")
+                dp_sel_ndx_file = op_in["task_path"].joinpath(dp_sel_ndx)
+                walker_idx = int(op_in["task_path"].name)
+                slice_dump(dump="out.dump",walker_idx = walker_idx,selected_idx=high_devi, output=sel_lmp_name,style="dpdata")
+                for ii, sel in enumerate(high_devi):
+                    dp_conf_list.append(op_in["task_path"].joinpath(sel_lmp_name.format(walker=walker_idx,idx=sel)))
+                
+                x_list = np.array([i for i in range(len(f_list))])
+                nframes = len(x_list)
+                groups = ["other"]*nframes
+                y_list = [int(j) for j in high_devi]
+                for j in y_list:
+                    groups[j] = "choice"
+                cdict = {"other": 'black', "choice": 'red'}
+                groups = np.array(groups)
+                fig, ax = plt.subplots()
+                for g in ["other","choice"]:
+                    ix = np.where(groups == g)
+                    ax.scatter(x_list[ix], f_list[ix], c = cdict[g], label = g, s = 20)
+                plt.title("max model devi of DP model")
+                plt.legend()
+                plt.xlabel("simu frames")
+                plt.ylabel("max model devi")
+                plt.savefig(dp_model_devi_fig)
+                dp_model_devi_file = op_in["task_path"].joinpath(dp_model_devi_file_name)
+                dp_model_devi_fig_file = op_in["task_path"].joinpath(dp_model_devi_fig)
+            
+            projected_fig_file = None
+            if "proj_info" in op_in["exploration_config"]:
+                proj_info = op_in["exploration_config"]["proj_info"]
+                # plot projected cv fig
+                proj_mode = proj_info["proj_mode"]
+                if proj_mode == "path":
+                    proj_cv_index = proj_info["proj_cv_index"]
+                    proj_lm = proj_info["path_lm"]
+                    path_list = proj_info["path_list"]
+                    cv_data = np.loadtxt("plm.out")[:,2:]
+                    proj_cv_data = cv_data[:,proj_cv_index]                                      
+                    slist = s(proj_cv_data, proj_lm, path_list)
+                    zlist = z(proj_cv_data, proj_lm, path_list)
+                    point_numbers = list(range(len(slist)))
+                    plt.figure(figsize=(10, 8), dpi=100)
+                    plt.scatter(slist,zlist, c = point_numbers)
+                    plt.xlabel("s")
+                    plt.ylabel("z")
+                    plt.title("Trajectory projection")
+                    plt.colorbar(label="steps")
+                    plt.savefig("path.png")
+                    projected_fig_file = op_in["task_path"].joinpath("path.png")
+            
             if op_in["models"] is not None:
                 # plot bias during simulation
                 bias = np.loadtxt(plumed_output_name)[:,1]
@@ -234,6 +320,11 @@ def execute(
                 "plm_out": op_in["task_path"].joinpath(plumed_output_name),
                 "bias_fig": bias_fig_file,
                 "model_devi_fig": model_devi_fig_file,
+                "dp_model_devi_fig": dp_model_devi_fig_file,
+                "dp_model_devi": dp_model_devi_file,
+                "dp_selected_indices": dp_sel_ndx_file,
+                "dp_selected_confs": dp_conf_list,
+                "projected_fig": projected_fig_file,
                 "md_log": op_in["task_path"].joinpath(mdrun_log),
                 "trajectory": op_in["task_path"].joinpath(traj_name),
                 "conf_out": op_in["task_path"].joinpath(conf_out),
diff --git a/rid/op/run_label.py b/rid/op/run_label.py
index a1edc951..0875cff5 100644
--- a/rid/op/run_label.py
+++ b/rid/op/run_label.py
@@ -7,7 +7,8 @@
     OPIO,
     OPIOSign,
     Artifact,
-    Parameter
+    Parameter,
+    BigParameter
 )
 from rid.constants import (
         gmx_conf_name,
@@ -54,10 +55,10 @@ def get_input_sign(cls):
             {
                 "task_path": Artifact(Path),
                 "forcefield": Artifact(Path, optional=True),
-                "label_config": Dict,
-                "cv_config": Dict,
+                "label_config": BigParameter(Dict),
+                "cv_config": BigParameter(Dict),
                 "at": Artifact(Path, optional=True),
-                "task_name": str,
+                "task_name": BigParameter(str),
                 "tail": Parameter(float, default=0.9),
                 "index_file": Artifact(Path, optional=True),
                 "dp_files": Artifact(List[Path], optional=True),
@@ -70,11 +71,11 @@ def get_input_sign(cls):
     def get_output_sign(cls):
         return OPIOSign(
             {
-                "plm_out": Artifact(Path),
-                "cv_forces": Artifact(Path),
-                "mf_info": Artifact(Path,optional=True),
-                "mf_fig": Artifact(Path),
-                "md_log": Artifact(Path)
+                "plm_out": Artifact(Path, archive = None),
+                "cv_forces": Artifact(Path, archive = None),
+                "mf_info": Artifact(Path,optional=True, archive = None),
+                "mf_fig": Artifact(Path, archive = None),
+                "md_log": Artifact(Path, archive = None)
             }
         )
 
diff --git a/rid/op/run_model_devi.py b/rid/op/run_model_devi.py
new file mode 100644
index 00000000..b7922865
--- /dev/null
+++ b/rid/op/run_model_devi.py
@@ -0,0 +1,128 @@
+from typing import List, Optional, Dict
+from pathlib import Path
+import numpy as np
+from dflow.python import (
+    OP,
+    OPIO,
+    OPIOSign,
+    Artifact,
+    Parameter,
+    BigParameter
+)
+from rid.utils import save_txt, set_directory
+from rid.constants import sel_gro_name, cv_init_label, model_devi_name, model_devi_precision, sel_ndx_name
+from rid.select.conf_select import select_from_devi
+from rid.common.mol import get_distance_from_atomid
+from rid.select.model_devi import make_std
+import os
+from matplotlib import pyplot as plt
+
+class RunModelDevi(OP):
+
+    """
+    `RunModelDevi` calculates model deviations for each chosen representive cluster frames from `PrepSelect` and select 
+    ones with high uncertainty from them.
+    As rid-kit is based on `Gromacs`, please provide trajectories in `.xtc` format (single-point precision) and NN models in 
+    `.pb` format. 
+    Warning: We highly recommend use `slice_mode = "gmx"` due to the inconsistent format convention of `mdtraj` that may lead to topology
+    mismatch of next label steps. If you use `mdtraj` mode, please make sure the name conventions of molecules in Gromacs topology 
+    files satisfy PDB standards. This could happen in some old Gromcas version.
+    """
+
+    @classmethod
+    def get_input_sign(cls):
+        return OPIOSign(
+            {
+                "models": Artifact(List[Path]),
+                "plm_out": Artifact(Path),
+                "selected_indices": Artifact(Path, optional=True),
+                "trust_lvl_1": float,
+                "task_name": str,
+                "block_tag": str
+            }
+        )
+
+    @classmethod
+    def get_output_sign(cls):
+        return OPIOSign(
+            {
+                "model_devi": Artifact(Path, optional=True, archive = None),
+                "model_devi_png": Artifact(Path, optional=True, archive = None)
+            }
+        )
+
+    @OP.exec_sign_check
+    def execute(
+        self,
+        op_in: OPIO,
+    ) -> OPIO:
+
+        r"""Execute the OP.
+        
+        Parameters
+        ----------
+        op_in : dict
+            Input dict with components:
+
+            - `models`: Artifact(List[Path], optional=True),
+            - `plm_out`: Artifact(Path)
+          
+        Returns
+        -------
+            Output dict with components:
+        
+            - `model_devi`: (`Artifact(Path)`) Model deviation values of selected conformation files (`selected_confs`) .
+        """
+        task_path = Path(op_in["task_name"])
+        task_path.mkdir(exist_ok=True, parents=True)
+
+        with set_directory(task_path):
+            # if "units" in cv_config:
+            #     length_units = cv_config["units"]
+            # else:
+            #     length_units = None
+            # if length_units == None or length_units == "nm":
+            #     cv_data = np.loadtxt(op_in["plm_out"])[:,2:]
+            # elif length_units == "A" or length_units == "Angstrom":
+            #     cv_data = 10*np.loadtxt(op_in["plm_out"])[:,2:]
+            # else:
+            #     raise ValueError("Valid length units must be nm or A")
+            trust_lvl = op_in["trust_lvl_1"]
+            cv_data = np.loadtxt(op_in["plm_out"])[:,2:]
+            stds = make_std(cv_data, models=op_in["models"])
+            save_txt("cls_"+model_devi_name, stds, fmt=model_devi_precision)
+            
+            _selected_idx = np.load(op_in["selected_indices"])
+     
+            nframes = cv_data.shape[0]
+            groups = ["other"]*nframes
+            xlist = np.linspace(0,nframes-1,nframes,dtype=int)
+            ylist = [int(j) for j in _selected_idx]                
+            for y_ in ylist:
+                groups[y_] = "train"
+            cdict = {"other": 'black',"train": 'red'}
+            fig, ax = plt.subplots()
+            groups = np.array(groups)
+            
+            for g in ["other","train"]:
+                ix = np.where(groups == g)
+                ax.scatter(xlist[ix], stds[ix], c = cdict[g], label = g, s = 20)
+            plt.axhline(y=trust_lvl, color='b', linestyle='-',label="trust level 1")
+            plt.xlabel("configuration list")
+            plt.ylabel("force std (KJ/(mol*length))")
+            plt.legend()
+            plt.title("prediction error")
+            fig.set_size_inches(8, 6, forward=True)
+            plt.savefig("model_devi.png")
+                
+        model_devi_fig_path = None
+        if os.path.exists(task_path.joinpath("model_devi.png")):
+            model_devi_fig_path = task_path.joinpath("model_devi.png")
+            
+        op_out = OPIO(
+            {
+               "model_devi": task_path.joinpath("cls_"+model_devi_name),
+               "model_devi_png": model_devi_fig_path
+            }
+        )
+        return op_out
\ No newline at end of file
diff --git a/rid/op/run_select.py b/rid/op/run_select.py
index 79eacf8b..24ab1e91 100644
--- a/rid/op/run_select.py
+++ b/rid/op/run_select.py
@@ -6,13 +6,16 @@
     OPIO,
     OPIOSign,
     Artifact,
-    Parameter
+    Parameter,
+    BigParameter
 )
 from rid.utils import save_txt, set_directory
 from rid.constants import sel_gro_name, sel_lmp_name, cv_init_label, model_devi_name, model_devi_precision, sel_ndx_name
 from rid.select.conf_select import select_from_devi
-from rid.common.mol import slice_xtc, slice_dump
+from rid.common.mol import slice_xtc
+from rid.common.mol_dpdata import slice_dump
 from rid.select.model_devi import make_std
+import json
 
 
 class RunSelect(OP):
@@ -39,6 +42,8 @@ def get_input_sign(cls):
                 "trust_lvl_2": float,
                 "xtc_traj": Artifact(Path),
                 "topology": Artifact(Path),
+                "label_config": BigParameter(Dict),
+                "type_map": Parameter(Optional[List], default=[]),
                 "dt": Parameter(Optional[float], default=None),
                 "output_freq": Parameter(Optional[float], default=None),
                 "slice_mode": Parameter(str, default="gmx")
@@ -49,11 +54,11 @@ def get_input_sign(cls):
     def get_output_sign(cls):
         return OPIOSign(
             {
-                "selected_confs": Artifact(List[Path]),
-                "selected_cv_init": Artifact(List[Path]),
-                "model_devi": Artifact(Path, optional=True),
-                "selected_indices": Artifact(Path),
-                "selected_conf_tags": Dict
+                "selected_confs": Artifact(List[Path], archive = None),
+                "selected_cv_init": Artifact(List[Path], archive = None),
+                "model_devi": Artifact(Path, optional=True, archive = None),
+                "selected_indices": Artifact(Path, archive = None),
+                "selected_conf_tags": Artifact(Path, archive= None)
             }
         )
 
@@ -107,7 +112,7 @@ def execute(
                 save_txt("cls_"+model_devi_name, stds, fmt=model_devi_precision)
                 _selected_idx = select_from_devi(stds, op_in["trust_lvl_1"])
             sel_idx = cls_sel_idx[_selected_idx]
-            save_txt(sel_ndx_name, sel_idx, fmt="%d")
+            np.save(sel_ndx_name, sel_idx)
             sel_data = cls_sel_data[_selected_idx]
             if op_in["slice_mode"] == "gmx":
                 assert op_in["dt"] is not None, "Please provide time step to slice trajectory."
@@ -119,7 +124,12 @@ def execute(
                 slice_xtc(xtc=op_in["xtc_traj"], top=op_in["topology"],
                         walker_idx = walker_idx, selected_idx=sel_idx, output=sel_gro_name, style="mdtraj")
             elif op_in["slice_mode"] == "dpdata":
-                slice_dump(dump=op_in["xtc_traj"],walker_idx = walker_idx,selected_idx=sel_idx, output=sel_lmp_name, style="dpdata")
+                if op_in["label_config"]["type"] == "gmx":
+                    slice_dump(dump=op_in["xtc_traj"],walker_idx = walker_idx,selected_idx=sel_idx, output=sel_gro_name, type_map = op_in["type_map"], style="dpdata")
+                elif op_in["label_config"]["type"] == "lmp":
+                    slice_dump(dump=op_in["xtc_traj"],walker_idx = walker_idx,selected_idx=sel_idx, output=sel_lmp_name, type_map = op_in["type_map"],style="dpdata")
+                else:
+                    raise ValueError("Invalid labeling type, only support gmx and lmp")
             else:
                 raise RuntimeError("Unknown Style for Slicing Trajectory.")
             conf_list = []
@@ -127,13 +137,19 @@ def execute(
             conf_tags = {}
             for ii, sel in enumerate(sel_idx):
                 if op_in["slice_mode"] == "dpdata":
-                    conf_list.append(task_path.joinpath(sel_lmp_name.format(walker=walker_idx,idx=sel)))
-                    conf_tags[sel_lmp_name.format(walker = walker_idx,idx=sel)] = f"{op_in['task_name']}_{sel}"
+                    if op_in["label_config"]["type"] == "gmx":
+                        conf_list.append(task_path.joinpath(sel_gro_name.format(walker=walker_idx,idx=sel)))
+                        conf_tags[sel_gro_name.format(walker = walker_idx,idx=sel)] = f"{op_in['task_name']}_{sel}"
+                    elif op_in["label_config"]["type"] == "lmp":
+                        conf_list.append(task_path.joinpath(sel_lmp_name.format(walker=walker_idx,idx=sel)))
+                        conf_tags[sel_lmp_name.format(walker = walker_idx,idx=sel)] = f"{op_in['task_name']}_{sel}"
                 elif op_in["slice_mode"] == "gmx" or op_in["slice_mode"] == "mdtraj" :
                     conf_list.append(task_path.joinpath(sel_gro_name.format(walker=walker_idx,idx=sel)))
                     conf_tags[sel_gro_name.format(walker = walker_idx,idx=sel)] = f"{op_in['task_name']}_{sel}"
                 save_txt(cv_init_label.format(walker=walker_idx,idx=sel), sel_data[ii])
                 cv_init_list.append(task_path.joinpath(cv_init_label.format(walker=walker_idx,idx=sel)))
+            with open("conf.json", "w") as f:
+                json.dump(conf_tags,f)
             
         op_out = OPIO(
             {
@@ -141,7 +157,7 @@ def execute(
                "selected_cv_init": cv_init_list,
                "model_devi": task_path.joinpath("cls_"+model_devi_name),
                "selected_indices": task_path.joinpath(sel_ndx_name),
-               "selected_conf_tags": conf_tags
+               "selected_conf_tags": task_path.joinpath("conf.json")
             }
         )
         return op_out
\ No newline at end of file
diff --git a/rid/op/run_train.py b/rid/op/run_train.py
index 0e519150..39140a26 100644
--- a/rid/op/run_train.py
+++ b/rid/op/run_train.py
@@ -6,7 +6,8 @@
     OPIO,
     OPIOSign,
     Artifact,
-    Parameter
+    Parameter,
+    BigParameter
 )
 from rid.constants import tf_model_name, train_fig, train_log
 from rid.nn.train_net import train
@@ -28,7 +29,7 @@ def get_input_sign(cls):
                 "model_tag": str,
                 "angular_mask": List,
                 "data": Artifact(Path),
-                "train_config": Dict
+                "train_config": BigParameter(Dict)
             }
         )
 
diff --git a/rid/select/cluster.py b/rid/select/cluster.py
index e43924c0..c3119e67 100644
--- a/rid/select/cluster.py
+++ b/rid/select/cluster.py
@@ -105,6 +105,26 @@ def mk_cluster(dist, distance_threshold):
     return cluster.labels_
 
 
+def chooseClusterCenter(dist:np.ndarray, conf_ids:list):
+    id_min = conf_ids[0]
+    min_loss = None
+    # n_tot_frame = dist.shape[0]
+    for id in conf_ids:
+        rmsd = dist[id]
+        loss = 0
+        for i, conf in enumerate(conf_ids):
+            rms = rmsd[conf]
+            loss += rms * rms
+        loss = np.sqrt(loss)
+        if min_loss == None:
+            min_loss = loss
+        else:
+            if loss < min_loss:
+                id_min = id
+                min_loss = loss
+    return [id_min]
+
+
 def sel_from_cluster(cvs, threshold, task_path, angular_mask=None, weights=None, max_selection=1000):
     if len(cvs) <= 1:
         return cvs
@@ -133,7 +153,8 @@ def sel_from_cluster(cvs, threshold, task_path, angular_mask=None, weights=None,
     cls_sel = []
     np.random.seed(seed=None)
     for cluster, _ in cls_map:
-        _ret = np.random.choice(cluster, 1)
+        # _ret = np.random.choice(cluster, 1)
+        _ret = chooseClusterCenter(dist, cluster)
         cls_sel.append(_ret[0])
     if len(cls_sel) > max_selection:
         cls_sel = cls_sel[:max_selection]
diff --git a/rid/superop/blocks.py b/rid/superop/blocks.py
index 51ec9bac..a2e05fc8 100755
--- a/rid/superop/blocks.py
+++ b/rid/superop/blocks.py
@@ -40,7 +40,9 @@ def __init__(
         label_op: OP,
         data_op: OP,
         train_op: OP,
+        model_devi_op: OP,
         train_config: Dict,
+        model_devi_config: Dict,
         upload_python_package = None,
         retry_times = None
     ):
@@ -63,6 +65,7 @@ def __init__(
             "output_freq": InputParameter(type=float, value=2500),
             "slice_mode": InputParameter(type=str, value="gmx"),
             "label_config": InputParameter(type=Dict),
+            "type_map": InputParameter(type=List, value = []),
             "tail": InputParameter(type=float, value=0.9),
             "train_config": InputParameter(type=Dict)
         }        
@@ -107,7 +110,9 @@ def __init__(
             label_op,
             data_op,
             train_op,
+            model_devi_op,
             train_config,
+            model_devi_config,
             upload_python_package = upload_python_package,
             retry_times = retry_times
         )            
@@ -140,7 +145,9 @@ def _first_run_block(
         label_op: OP,
         data_op: OP,
         train_op: OP,
+        model_devi_op: OP,
         train_config : Dict,
+        model_devi_config: Dict,
         upload_python_package : str = None,
         retry_times: int = None
     ):
@@ -173,6 +180,7 @@ def _first_run_block(
         "Selection",
         template=select_op,
         parameters={
+            "label_config": block_steps.inputs.parameters["label_config"],
             "trust_lvl_1" : block_steps.inputs.parameters["trust_lvl_1"],
             "trust_lvl_2": block_steps.inputs.parameters["trust_lvl_2"],
             "cluster_threshold": block_steps.inputs.parameters["cluster_threshold"],
@@ -184,6 +192,7 @@ def _first_run_block(
             "dt": block_steps.inputs.parameters["dt"],
             "output_freq": block_steps.inputs.parameters["output_freq"],
             "slice_mode": block_steps.inputs.parameters["slice_mode"],
+            "type_map": block_steps.inputs.parameters["type_map"],
             "if_make_threshold": True,
             "task_names" : block_steps.inputs.parameters['walker_tags'],
             "block_tag" : block_steps.inputs.parameters['block_tag'],
@@ -205,7 +214,6 @@ def _first_run_block(
             "label_config": block_steps.inputs.parameters['label_config'],
             "cv_config": block_steps.inputs.parameters['cv_config'],
             "tail": block_steps.inputs.parameters['tail'],
-            "conf_tags" : selection.outputs.parameters['selected_conf_tags'],
             "block_tag" : block_steps.inputs.parameters['block_tag'],
         },
         artifacts={
@@ -216,7 +224,8 @@ def _first_run_block(
             "index_file": block_steps.inputs.artifacts['index_file'],
             "inputfile": block_steps.inputs.artifacts['inputfile'],
             "dp_files": block_steps.inputs.artifacts['dp_files'],
-            "cv_file": block_steps.inputs.artifacts['cv_file']
+            "cv_file": block_steps.inputs.artifacts['cv_file'],
+            "conf_tags" : selection.outputs.artifacts['selected_conf_tags']
         },
         key = '{}-label'.format(block_steps.inputs.parameters['block_tag'])
     )
@@ -262,6 +271,38 @@ def _first_run_block(
         **train_config,
     )
     block_steps.add(train)
+    
+    model_devi_config = deepcopy(model_devi_config)
+    model_devi_template_config = model_devi_config.pop('template_config')
+    model_devi_executor = init_executor(model_devi_config.pop('executor'))
+    deviation = Step(
+        "ModelDeviation",
+        template=PythonOPTemplate(
+            model_devi_op,
+            python_packages = upload_python_package,
+            retry_on_transient_error = retry_times,
+            slices=Slices("{{item}}",
+                input_parameter=["trust_lvl_1","task_name"],
+                input_artifact=["plm_out","selected_indices"],
+                output_artifact=["model_devi","model_devi_png"]),
+            **model_devi_template_config,
+        ),
+        parameters={
+            "trust_lvl_1": block_steps.inputs.parameters["trust_lvl_1"],
+            "block_tag" : block_steps.inputs.parameters['block_tag'],
+            "task_name": block_steps.inputs.parameters['walker_tags'],
+        },
+        artifacts={
+            "models" : train.outputs.artifacts["model"],
+            "plm_out": exploration.outputs.artifacts["plm_out"],
+            "selected_indices": selection.outputs.artifacts["selected_indices"]
+        },
+        executor = model_devi_executor,
+        with_param=argo_range(argo_len(block_steps.inputs.parameters["walker_tags"])),
+        key = '{}-model-devi'.format(block_steps.inputs.parameters['block_tag'])+"-{{item}}",
+        **model_devi_config
+    )
+    block_steps.add(deviation)
 
     block_steps.outputs.artifacts["models"]._from = train.outputs.artifacts["model"]
     block_steps.outputs.artifacts["data"]._from = gen_data.outputs.artifacts["data"]
@@ -288,8 +329,10 @@ def __init__(
         data_op: OP,
         adjust_lvl_op: OP,
         train_op: OP,  
+        model_devi_op: OP,
         adjust_lvl_config: Dict,
         train_config: Dict,
+        model_devi_config: Dict,
         upload_python_package = None,
         retry_times = None
     ):
@@ -315,6 +358,7 @@ def __init__(
             "label_config": InputParameter(type=Dict),
             "tail": InputParameter(type=float, value=0.9),
             "train_config": InputParameter(type=Dict),
+            "type_map": InputParameter(type=List, value=[]),
             "adjust_amplifier": InputParameter(type=float, value=1.5),
             "max_level_multiple": InputParameter(type=float, value=8.0),
         }        
@@ -363,8 +407,10 @@ def __init__(
             data_op,
             adjust_lvl_op,
             train_op,
+            model_devi_op,
             adjust_lvl_config,
             train_config,
+            model_devi_config,
             upload_python_package = upload_python_package,
             retry_times = retry_times
         )            
@@ -398,8 +444,10 @@ def _iter_block(
         data_op: OP,
         adjust_lvl_op: OP,
         train_op: OP,
+        model_devi_op: OP,
         adjust_lvl_config : Dict,
         train_config : Dict,
+        model_devi_config: Dict,
         upload_python_package : str = None,
         retry_times: int = None
     ):
@@ -433,6 +481,7 @@ def _iter_block(
         "Selection",
         template=select_op,
         parameters={
+            "label_config": block_steps.inputs.parameters["label_config"],
             "trust_lvl_1" : block_steps.inputs.parameters["init_trust_lvl_1"],
             "trust_lvl_2": block_steps.inputs.parameters["init_trust_lvl_2"],
             "cluster_threshold": block_steps.inputs.parameters["cluster_threshold"],
@@ -442,6 +491,7 @@ def _iter_block(
             "dt": block_steps.inputs.parameters["dt"],
             "output_freq": block_steps.inputs.parameters["output_freq"],
             "slice_mode": block_steps.inputs.parameters["slice_mode"],
+            "type_map": block_steps.inputs.parameters["type_map"],
             "if_make_threshold": False,
             "task_names" : block_steps.inputs.parameters['walker_tags'],
             "block_tag" : block_steps.inputs.parameters['block_tag'],
@@ -497,7 +547,6 @@ def _iter_block(
             "label_config": block_steps.inputs.parameters['label_config'],
             "cv_config": block_steps.inputs.parameters['cv_config'],
             "tail": block_steps.inputs.parameters['tail'],
-            "conf_tags" : selection.outputs.parameters['selected_conf_tags'],
             "block_tag" : block_steps.inputs.parameters['block_tag'],
         },
         artifacts={
@@ -508,7 +557,8 @@ def _iter_block(
             "at": selection.outputs.artifacts["selected_cv_init"],
             "index_file": block_steps.inputs.artifacts['index_file'],
             "dp_files": block_steps.inputs.artifacts['dp_files'],
-            "cv_file": block_steps.inputs.artifacts['cv_file']
+            "cv_file": block_steps.inputs.artifacts['cv_file'],
+            "conf_tags" : selection.outputs.artifacts['selected_conf_tags']
         },
         key = '{}-label'.format(block_steps.inputs.parameters['block_tag'])
     )
@@ -554,6 +604,38 @@ def _iter_block(
         **train_config,
     )
     block_steps.add(train)
+    
+    model_devi_config = deepcopy(model_devi_config)
+    model_devi_template_config = model_devi_config.pop('template_config')
+    model_devi_executor = init_executor(model_devi_config.pop('executor'))
+    deviation = Step(
+        "ModelDeviation",
+        template=PythonOPTemplate(
+            model_devi_op,
+            python_packages = upload_python_package,
+            retry_on_transient_error = retry_times,
+            slices=Slices("{{item}}",
+                input_parameter=["trust_lvl_1","task_name"],
+                input_artifact=["plm_out","selected_indices"],
+                output_artifact=["model_devi","model_devi_png"]),
+            **model_devi_template_config,
+        ),
+        parameters={
+            "trust_lvl_1": block_steps.inputs.parameters["trust_lvl_1"],
+            "block_tag" : block_steps.inputs.parameters['block_tag'],
+            "task_name": block_steps.inputs.parameters['walker_tags'],
+        },
+        artifacts={
+            "models" : train.outputs.artifacts["model"],
+            "plm_out": exploration.outputs.artifacts["plm_out"],
+            "selected_indices": selection.outputs.artifacts["selected_indices"]
+        },
+        executor = model_devi_executor,
+        with_param=argo_range(argo_len(block_steps.inputs.parameters["walker_tags"])),
+        key = '{}-model-devi'.format(block_steps.inputs.parameters['block_tag'])+"-{{item}}",
+        **model_devi_config
+    )
+    block_steps.add(deviation)
 
     block_steps.outputs.artifacts["models"]._from = train.outputs.artifacts["model"]
     block_steps.outputs.artifacts["data"]._from = gen_data.outputs.artifacts["data"]
diff --git a/rid/superop/exploration.py b/rid/superop/exploration.py
index 7a85eb5c..7d45071b 100755
--- a/rid/superop/exploration.py
+++ b/rid/superop/exploration.py
@@ -128,7 +128,11 @@ def _exploration(
     run_template_config = run_config.pop('template_config')
     prep_executor = init_executor(prep_config.pop('executor'))
     run_executor = init_executor(run_config.pop('executor'))
-    prep_exploration = Step(
+    prep_merge = False
+    if prep_executor is not None:
+        prep_merge = prep_executor.merge_sliced_step
+    if prep_merge:
+        prep_exploration = Step(
         'prep-exploration',
         template=PythonOPTemplate(
             prep_exploration_op,
@@ -159,11 +163,48 @@ def _exploration(
         key = step_keys["prep_exploration"]+"-{{item}}",
         with_param=argo_range(argo_len(exploration_steps.inputs.parameters['task_names'])),
         executor = prep_executor,
-        **prep_config,
+        **prep_config
+    )
+    else:
+        prep_exploration = Step(
+            'prep-exploration',
+            template=PythonOPTemplate(
+                prep_exploration_op,
+                python_packages = upload_python_package,
+                retry_on_transient_error = retry_times,
+                slices=Slices(sub_path = True,
+                    input_parameter=["task_name", "trust_lvl_1", "trust_lvl_2"],
+                    input_artifact=["conf"],
+                    output_artifact=["task_path"],
+                    output_parameter=["cv_dim"]
+                ),
+                **prep_template_config,
+            ),
+            parameters={
+                "trust_lvl_1" : exploration_steps.inputs.parameters['trust_lvl_1'],
+                "trust_lvl_2": exploration_steps.inputs.parameters['trust_lvl_2'],
+                "exploration_config" : exploration_steps.inputs.parameters['exploration_config'],
+                "cv_config" : exploration_steps.inputs.parameters['cv_config'],
+                "task_name": exploration_steps.inputs.parameters['task_names'],
+                "block_tag": exploration_steps.inputs.parameters["block_tag"]
+            },
+            artifacts={
+                "models" : exploration_steps.inputs.artifacts['models'],
+                "topology" :exploration_steps.inputs.artifacts['topology'],
+                "conf" : exploration_steps.inputs.artifacts['confs'],
+                "cv_file": exploration_steps.inputs.artifacts['cv_file']
+            },
+            key = step_keys["prep_exploration"]+"-{{item.order}}",
+            executor = prep_executor,
+            **prep_config
     )
     exploration_steps.add(prep_exploration)
 
-    run_exploration = Step(
+    run_merge = False
+    if run_executor is not None:
+        run_merge = run_executor.merge_sliced_step
+    if run_merge:
+        run_exploration = Step(
         'run-exploration',
         template=PythonOPTemplate(
             run_exploration_op,
@@ -171,7 +212,7 @@ def _exploration(
             retry_on_transient_error = retry_times,
             slices=Slices("{{item}}",
                 input_artifact=["task_path"],
-                output_artifact=["plm_out", "bias_fig","model_devi_fig", "trajectory", "md_log", "conf_out"]
+                output_artifact=["plm_out", "bias_fig","model_devi_fig", "dp_model_devi_fig", "dp_model_devi", "dp_selected_indices","dp_selected_confs","projected_fig","trajectory", "md_log", "conf_out"]
             ),
             **run_template_config,
         ),
@@ -190,8 +231,37 @@ def _exploration(
         key = step_keys["run_exploration"]+"-{{item}}",
         executor = run_executor,
         with_param=argo_range(argo_len(exploration_steps.inputs.parameters['task_names'])),
-        **run_config,
+        **run_config
     )
+    else:
+        run_exploration = Step(
+            'run-exploration',
+            template=PythonOPTemplate(
+                run_exploration_op,
+                python_packages = upload_python_package,
+                retry_on_transient_error = retry_times,
+                slices=Slices(sub_path = True,
+                    input_artifact=["task_path"],
+                    output_artifact=["plm_out", "bias_fig","model_devi_fig","dp_model_devi_fig", "dp_model_devi", "dp_selected_indices","dp_selected_confs","projected_fig","trajectory", "md_log", "conf_out"]
+                ),
+                **run_template_config,
+            ),
+            parameters={
+                "exploration_config" : exploration_steps.inputs.parameters["exploration_config"]
+            },
+            artifacts={
+                "task_path" : prep_exploration.outputs.artifacts["task_path"],
+                "forcefield": exploration_steps.inputs.artifacts['forcefield'],
+                "models" : exploration_steps.inputs.artifacts['models'],
+                "index_file": exploration_steps.inputs.artifacts['index_file'],
+                "dp_files": exploration_steps.inputs.artifacts['dp_files'],
+                "cv_file": exploration_steps.inputs.artifacts['cv_file'],
+                "inputfile": exploration_steps.inputs.artifacts['inputfile']
+            },
+            key = step_keys["run_exploration"]+"-{{item.order}}",
+            executor = run_executor,
+            **run_config,
+        )
     exploration_steps.add(run_exploration)
 
     exploration_steps.outputs.parameters["cv_dim"].value_from_parameter = prep_exploration.outputs.parameters["cv_dim"]
diff --git a/rid/superop/label.py b/rid/superop/label.py
index 48545ca8..e424526e 100755
--- a/rid/superop/label.py
+++ b/rid/superop/label.py
@@ -40,7 +40,6 @@ def __init__(
             "label_config": InputParameter(type=Dict),
             "cv_config": InputParameter(type=Dict),
             "tail": InputParameter(type=float, value=0.9),
-            "conf_tags": InputParameter(type=List),
             "block_tag" : InputParameter(type=str, value="")
         }        
         self._input_artifacts = {
@@ -52,7 +51,8 @@ def __init__(
             "at": InputArtifact(optional=True),
             "index_file": InputArtifact(optional=True),
             "dp_files": InputArtifact(optional=True),
-            "cv_file": InputArtifact(optional=True)
+            "cv_file": InputArtifact(optional=True),
+            "conf_tags": InputArtifact(optional=True)
         }
         self._output_parameters = {
         }
@@ -140,10 +140,9 @@ def _label(
             retry_on_transient_error = retry_times,
             **prep_template_config,
         ),
-        parameters={
-            "conf_tags": label_steps.inputs.parameters['conf_tags'],  
-        },
+        parameters={},
         artifacts={
+            "conf_tags": label_steps.inputs.artifacts['conf_tags'],  
             "confs": label_steps.inputs.artifacts['confs'],
         },
         key = step_keys['check_label_inputs'],
@@ -152,7 +151,11 @@ def _label(
     )
     label_steps.add(check_label_inputs)
 
-    prep_label = Step(
+    prep_merge = False
+    if prep_executor is not None:
+        prep_merge = prep_executor.merge_sliced_step
+    if prep_merge:
+        prep_label = Step(
         'prep-label',
         template=PythonOPTemplate(
             prep_label_op,
@@ -179,17 +182,89 @@ def _label(
         executor = prep_executor,
         with_param=argo_range(argo_len(check_label_inputs.outputs.parameters['conf_tags'])),
         when = "%s > 0" % (check_label_inputs.outputs.parameters["if_continue"]),
-        **prep_config,
+        **prep_config
+    )
+    else:
+        prep_label = Step(
+        'prep-label',
+        template=PythonOPTemplate(
+            prep_label_op,
+            python_packages = upload_python_package,
+            retry_on_transient_error = retry_times,
+            slices=Slices("{{item}}",
+                group_size=10,
+                pool_size=1,
+                input_parameter=["task_name"],
+                input_artifact=["conf", "at"],
+                output_artifact=["task_path"]),
+            **prep_template_config,
+        ),
+        parameters={
+            "label_config": label_steps.inputs.parameters['label_config'],
+            "cv_config": label_steps.inputs.parameters['cv_config'],
+            "task_name": check_label_inputs.outputs.parameters['conf_tags']
+        },
+        artifacts={
+            "topology": label_steps.inputs.artifacts['topology'],
+            "conf": label_steps.inputs.artifacts['confs'],
+            "at": label_steps.inputs.artifacts['at'],
+            "cv_file": label_steps.inputs.artifacts['cv_file']
+        },
+        key = step_keys['prep_label']+"-{{item}}",
+        executor = prep_executor,
+        with_param=argo_range(argo_len(check_label_inputs.outputs.parameters['conf_tags'])),
+        when = "%s > 0" % (check_label_inputs.outputs.parameters["if_continue"]),
+        **prep_config
     )
     label_steps.add(prep_label)
 
-    run_label = Step(
+    run_merge = False
+    if run_executor is not None:
+        run_merge = run_executor.merge_sliced_step
+    if run_merge:
+        run_label = Step(
+        'run-label',
+        template=PythonOPTemplate(
+            run_label_op,
+            python_packages = upload_python_package,
+            retry_on_transient_error = retry_times,
+            slices=Slices("{{item}}",
+                input_parameter=["task_name"],
+                input_artifact=["task_path","at"],
+                output_artifact=["plm_out","cv_forces","mf_info","mf_fig","md_log"]),
+            **run_template_config,
+        ),
+        parameters={
+            "label_config": label_steps.inputs.parameters["label_config"],
+            "cv_config": label_steps.inputs.parameters['cv_config'],
+            "task_name": check_label_inputs.outputs.parameters['conf_tags'],
+            "tail": label_steps.inputs.parameters['tail']
+        },
+        artifacts={
+            "forcefield": label_steps.inputs.artifacts['forcefield'],
+            "task_path": prep_label.outputs.artifacts["task_path"],
+            "index_file": label_steps.inputs.artifacts['index_file'],
+            "dp_files": label_steps.inputs.artifacts['dp_files'],
+            "cv_file": label_steps.inputs.artifacts['cv_file'],
+            "inputfile": label_steps.inputs.artifacts['inputfile'],
+            "at": label_steps.inputs.artifacts['at']
+        },
+        key = step_keys['run_label']+"-{{item}}",
+        executor = run_executor,
+        with_param=argo_range(argo_len(check_label_inputs.outputs.parameters['conf_tags'])),
+        continue_on_success_ratio = 0.75,
+        **run_config
+    )
+    else:
+        run_label = Step(
         'run-label',
         template=PythonOPTemplate(
             run_label_op,
             python_packages = upload_python_package,
             retry_on_transient_error = retry_times,
             slices=Slices("{{item}}",
+                group_size=10,
+                pool_size=1,
                 input_parameter=["task_name"],
                 input_artifact=["task_path","at"],
                 output_artifact=["plm_out","cv_forces","mf_info","mf_fig","md_log"]),
@@ -214,7 +289,7 @@ def _label(
         executor = run_executor,
         with_param=argo_range(argo_len(check_label_inputs.outputs.parameters['conf_tags'])),
         continue_on_success_ratio = 0.75,
-        **run_config,
+        **run_config
     )
     label_steps.add(run_label)
 
diff --git a/rid/superop/mcmc.py b/rid/superop/mcmc.py
new file mode 100644
index 00000000..2cdea066
--- /dev/null
+++ b/rid/superop/mcmc.py
@@ -0,0 +1,198 @@
+from typing import Dict, List
+from copy import deepcopy
+from dflow import (
+    InputParameter,
+    OutputParameter,
+    Inputs,
+    InputArtifact,
+    Outputs,
+    OutputArtifact,
+    Step,
+    Steps,
+    argo_range,
+    argo_len,
+)
+from dflow.python import(
+    PythonOPTemplate,
+    OP,
+    Slices
+)
+from rid.utils import init_executor
+
+
+class MCMC(Steps):
+    
+    r"""" MCMC SuperOP.
+    This SuperOP combines MCMC_Run OP and MCMC_Plot OP.
+    """
+    def __init__(
+        self,
+        name: str,
+        mcmc_run_op: OP,
+        mcmc_plot_op: OP,
+        run_config: Dict,
+        plot_config: Dict,
+        upload_python_package = None,
+        retry_times = None
+    ):
+        self._input_parameters = {
+            "mcmc_config" : InputParameter(type=Dict),
+            "task_names" : InputParameter(type=List[str]),
+            "block_tag" : InputParameter(type=str, value="")
+        }        
+        self._input_artifacts = {
+            "models" : InputArtifact()
+        }
+        self._output_parameters = {}
+        self._output_artifacts = {
+            "mcmc_fig": OutputArtifact()
+        }
+
+        super().__init__(        
+                name=name,
+                inputs=Inputs(
+                    parameters=self._input_parameters,
+                    artifacts=self._input_artifacts
+                ),
+                outputs=Outputs(
+                    parameters=self._output_parameters,
+                    artifacts=self._output_artifacts
+                ),
+            )
+
+        step_keys = {
+            "mcmc_run": "{}-mcmc-run".format(self.inputs.parameters["block_tag"]),
+            "mcmc_plot": "{}-mcmc-plot".format(self.inputs.parameters["block_tag"]),
+        }
+
+        self = _mcmc(
+            self, 
+            step_keys,
+            mcmc_run_op,
+            mcmc_plot_op,
+            run_config = run_config,
+            plot_config = plot_config,
+            upload_python_package = upload_python_package,
+            retry_times = retry_times
+        )            
+    
+    @property
+    def input_parameters(self):
+        return self._input_parameters
+
+    @property
+    def input_artifacts(self):
+        return self._input_artifacts
+
+    @property
+    def output_parameters(self):
+        return self._output_parameters
+
+    @property
+    def output_artifacts(self):
+        return self._output_artifacts
+
+    @property
+    def keys(self):
+        return self._keys
+
+
+def _mcmc(
+        mcmc_steps,
+        step_keys,
+        mcmc_run_op : OP,
+        mcmc_plot_op : OP,
+        run_config : Dict,
+        plot_config : Dict,
+        upload_python_package : str = None,
+        retry_times: int = None
+    ):
+    run_config = deepcopy(run_config)
+    plot_config = deepcopy(plot_config)
+    
+    run_template_config = run_config.pop('template_config')
+    plot_template_config = plot_config.pop('template_config')
+    run_executor = init_executor(run_config.pop('executor'))
+    plot_executor = init_executor(plot_config.pop('executor'))
+    
+    run_merge = False
+    if run_executor is not None:
+        run_merge = run_executor.merge_sliced_step
+    if run_merge:
+        mcmc_run = Step(
+        'mcmc-run',
+        template=PythonOPTemplate(
+            mcmc_run_op,
+            python_packages = upload_python_package,
+            retry_on_transient_error = retry_times,
+            slices=Slices("{{item}}",
+                input_parameter=["task_names"],
+                input_artifact=["models"],
+                output_artifact=["mcmc_1cv", "mcmc_2cv"]
+            ),
+            **run_template_config,
+        ),
+        parameters={
+            "mcmc_config" : mcmc_steps.inputs.parameters['mcmc_config'],
+            "task_names": mcmc_steps.inputs.parameters['task_names']
+        },
+        artifacts={
+            "models" : mcmc_steps.inputs.artifacts['models']
+        },
+        key = step_keys["mcmc_run"]+"-{{item}}",
+        with_param=argo_range(argo_len(mcmc_steps.inputs.parameters['task_names'])),
+        executor = run_executor,
+        **run_config
+    )
+    else:
+        mcmc_run = Step(
+            'mcmc-run',
+            template=PythonOPTemplate(
+                mcmc_run_op,
+                python_packages = upload_python_package,
+                retry_on_transient_error = retry_times,
+                slices=Slices(sub_path = True,
+                    input_parameter=["task_names"],
+                    input_artifact=["models"],
+                    output_artifact=["mcmc_1cv", "mcmc_2cv"]
+                ),
+                **run_template_config,
+            ),
+            parameters={
+                "mcmc_config" : mcmc_steps.inputs.parameters['mcmc_config'],
+                "task_names": mcmc_steps.inputs.parameters['task_names']
+            },
+            artifacts={
+                "models" : mcmc_steps.inputs.artifacts['models']
+            },
+            key = step_keys["mcmc_run"]+"-{{item.order}}",
+            executor = run_executor,
+            **run_config
+    )
+    mcmc_steps.add(mcmc_run)
+
+    mcmc_plot = Step(
+        'mcmc-plot',
+        template=PythonOPTemplate(
+            mcmc_plot_op,
+            python_packages = upload_python_package,
+            retry_on_transient_error = retry_times,
+            **plot_template_config,
+        ),
+        parameters={
+            "mcmc_config" : mcmc_steps.inputs.parameters["mcmc_config"]
+        },
+        artifacts={
+            "mcmc_1cv": mcmc_run.outputs.artifacts['mcmc_1cv'],
+            "mcmc_2cv": mcmc_run.outputs.artifacts['mcmc_2cv']
+        },
+        key = step_keys["mcmc_plot"],
+        executor = plot_executor,
+        **run_config
+    )
+
+    mcmc_steps.add(mcmc_plot)
+
+    mcmc_steps.outputs.artifacts["mcmc_fig"]._from = mcmc_plot.outputs.artifacts["mcmc_fig"]
+    
+    return mcmc_steps
\ No newline at end of file
diff --git a/rid/superop/selector.py b/rid/superop/selector.py
index 8b819b0d..8a77c4b9 100755
--- a/rid/superop/selector.py
+++ b/rid/superop/selector.py
@@ -37,6 +37,7 @@ def __init__(
         retry_times = None
     ):
         self._input_parameters = {
+            "label_config": InputParameter(type=Dict),
             "trust_lvl_1" : InputParameter(type=List[float], value=2.0),
             "trust_lvl_2": InputParameter(type=List[float], value=3.0),
             "cluster_threshold": InputParameter(type=List[float], value=1.0),
@@ -48,6 +49,7 @@ def __init__(
             "dt": InputParameter(type=float, value=0.02),
             "output_freq": InputParameter(type=float, value=2500),
             "slice_mode": InputParameter(type=str, value="gmx"),
+            "type_map": InputParameter(type=List, value=[]),
             "if_make_threshold": InputParameter(type=bool, value=False),
             "task_names" : InputParameter(type=List[str]),
             "block_tag" : InputParameter(type=str, value="")
@@ -60,15 +62,15 @@ def __init__(
         }
         self._output_parameters = {
             "cluster_threshold": OutputParameter(type=List[int]),
-            "numb_cluster": OutputParameter(type=List[int]),
-            "selected_conf_tags": OutputParameter(type=List),
+            "numb_cluster": OutputParameter(type=List[int])
         }
         self._output_artifacts = {
             "cluster_selection_index": OutputArtifact(),
             "selected_confs": OutputArtifact(),
             "selected_cv_init": OutputArtifact(),
             "model_devi": OutputArtifact(),
-            "selected_indices": OutputArtifact()
+            "selected_indices": OutputArtifact(),
+            "selected_conf_tags": OutputArtifact()
         }
 
         super().__init__(        
@@ -138,7 +140,11 @@ def _select(
     prep_executor = init_executor(prep_config.pop('executor'))
     run_executor = init_executor(run_config.pop('executor'))
 
-    prep_select = Step(
+    prep_merge = False
+    if prep_executor is not None:
+        prep_merge = prep_executor.merge_sliced_step
+    if prep_merge:
+        prep_select = Step(
         'prep-select',
         template=PythonOPTemplate(
             prep_select_op,
@@ -170,11 +176,48 @@ def _select(
         key = step_keys["prep_select"]+"-{{item}}",
         executor = prep_executor,
         with_param=argo_range(argo_len(select_steps.inputs.parameters['task_names'])),
+        **prep_config
+    )  
+    else:
+        prep_select = Step(
+        'prep-select',
+        template=PythonOPTemplate(
+            prep_select_op,
+            python_packages = upload_python_package,
+            retry_on_transient_error = retry_times,
+            slices=Slices(sub_path = True,
+                input_parameter=["cluster_threshold", "task_name"],
+                input_artifact=["plm_out"],
+                output_artifact=["cluster_fig","cluster_selection_index", "cluster_selection_data"],
+                output_parameter=["cluster_threshold", "numb_cluster"]
+                ),
+                
+            **prep_template_config,
+        ),
+        parameters={
+            "cluster_threshold": select_steps.inputs.parameters['cluster_threshold'],
+            "angular_mask": select_steps.inputs.parameters['angular_mask'],
+            "weights": select_steps.inputs.parameters['weights'],
+            "numb_cluster_upper": select_steps.inputs.parameters['numb_cluster_upper'],
+            "numb_cluster_lower": select_steps.inputs.parameters['numb_cluster_lower'],
+            "max_selection": select_steps.inputs.parameters['max_selection'],
+            "if_make_threshold": select_steps.inputs.parameters['if_make_threshold'],
+            "task_name": select_steps.inputs.parameters['task_names']
+        },
+        artifacts={
+            "plm_out": select_steps.inputs.artifacts['plm_out']
+        },
+        key = step_keys["prep_select"]+"-{{item.order}}",
+        executor = prep_executor,
         **prep_config,
     )
     select_steps.add(prep_select)
 
-    run_select = Step(
+    run_merge = False
+    if run_executor is not None:
+        run_merge = run_executor.merge_sliced_step
+    if run_merge:
+        run_select = Step(
         'run-select',
         template=PythonOPTemplate(
             run_select_op,
@@ -184,17 +227,18 @@ def _select(
                 "int({{item}})",
                 input_parameter=["task_name", "trust_lvl_1", "trust_lvl_2"],
                 input_artifact=["cluster_selection_index", "cluster_selection_data", "xtc_traj", "topology"],
-                output_artifact=["selected_confs", "selected_cv_init", "model_devi", "selected_indices"],
-                output_parameter=["selected_conf_tags"]
+                output_artifact=["selected_confs", "selected_cv_init", "model_devi", "selected_indices","selected_conf_tags"]
             ),
             **run_template_config,
         ),
         parameters={
+            "label_config":  select_steps.inputs.parameters["label_config"],
             "trust_lvl_1": select_steps.inputs.parameters["trust_lvl_1"],
             "trust_lvl_2": select_steps.inputs.parameters["trust_lvl_2"],
             "dt": select_steps.inputs.parameters["dt"],
             "output_freq": select_steps.inputs.parameters["output_freq"],
             "slice_mode": select_steps.inputs.parameters["slice_mode"],
+            "type_map": select_steps.inputs.parameters["type_map"],
             "task_name": select_steps.inputs.parameters['task_names']
         },
         artifacts={
@@ -207,14 +251,49 @@ def _select(
         key = step_keys["run_select"]+"-{{item}}",
         executor = run_executor,
         with_param=argo_range(argo_len(select_steps.inputs.parameters["task_names"])),
+        **run_config
+    )
+    else:
+        run_select = Step(
+        'run-select',
+        template=PythonOPTemplate(
+            run_select_op,
+            python_packages = upload_python_package,
+            retry_on_transient_error = retry_times,
+            slices=Slices(sub_path = True,
+                input_parameter=["task_name", "trust_lvl_1", "trust_lvl_2"],
+                input_artifact=["cluster_selection_index", "cluster_selection_data", "xtc_traj", "topology"],
+                output_artifact=["selected_confs", "selected_cv_init", "model_devi", "selected_indices","selected_conf_tags"]
+            ),
+            **run_template_config,
+        ),
+        parameters={
+            "label_config":  select_steps.inputs.parameters["label_config"],
+            "trust_lvl_1": select_steps.inputs.parameters["trust_lvl_1"],
+            "trust_lvl_2": select_steps.inputs.parameters["trust_lvl_2"],
+            "dt": select_steps.inputs.parameters["dt"],
+            "output_freq": select_steps.inputs.parameters["output_freq"],
+            "slice_mode": select_steps.inputs.parameters["slice_mode"],
+            "type_map": select_steps.inputs.parameters["type_map"],
+            "task_name": select_steps.inputs.parameters['task_names']
+        },
+        artifacts={
+            "cluster_selection_index": prep_select.outputs.artifacts["cluster_selection_index"],
+            "cluster_selection_data": prep_select.outputs.artifacts["cluster_selection_data"],
+            "models": select_steps.inputs.artifacts["models"],
+            "xtc_traj": select_steps.inputs.artifacts["xtc_traj"],
+            "topology": select_steps.inputs.artifacts["topology"]
+        },
+        key = step_keys["run_select"]+"-{{item.order}}",
+        executor = run_executor,
         **run_config,
     )
     select_steps.add(run_select)
 
     select_steps.outputs.parameters["cluster_threshold"].value_from_parameter = prep_select.outputs.parameters["cluster_threshold"]
     select_steps.outputs.parameters["numb_cluster"].value_from_parameter = prep_select.outputs.parameters["numb_cluster"]
-    select_steps.outputs.parameters["selected_conf_tags"].value_from_parameter = run_select.outputs.parameters["selected_conf_tags"]
 
+    select_steps.outputs.artifacts["selected_conf_tags"]._from = run_select.outputs.artifacts["selected_conf_tags"]
     select_steps.outputs.artifacts["cluster_selection_index"]._from = prep_select.outputs.artifacts["cluster_selection_index"]
     select_steps.outputs.artifacts["selected_confs"]._from = run_select.outputs.artifacts["selected_confs"]
     select_steps.outputs.artifacts["selected_cv_init"]._from = run_select.outputs.artifacts["selected_cv_init"]
diff --git a/rid/template/machine_bohrium.json b/rid/template/machine_bohrium.json
index 888069de..95e4e053 100755
--- a/rid/template/machine_bohrium.json
+++ b/rid/template/machine_bohrium.json
@@ -2,13 +2,13 @@
     "resources": {
         "bohrium1": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
             "machine_dict":{
-                "batch_type": "DpCloudServer",
-                "context_type": "DpCloudServerContext",
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
                 "local_root" : "./",
                 "remote_profile":{
                     "email": "",
@@ -32,13 +32,13 @@
         },
         "bohrium2": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
             "machine_dict":{
-                "batch_type": "DpCloudServer",
-                "context_type": "DpCloudServerContext",
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
                 "local_root" : "./",
                 "remote_profile":{
                     "email": "",
@@ -62,13 +62,13 @@
         },
         "bohrium3": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
             "machine_dict":{
-                "batch_type": "DpCloudServer",
-                "context_type": "DpCloudServerContext",
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
                 "local_root" : "./",
                 "remote_profile":{
                     "email": "",
@@ -92,13 +92,13 @@
         },
         "bohrium4": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
             "machine_dict":{
-                "batch_type": "DpCloudServer",
-                "context_type": "DpCloudServerContext",
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
                 "local_root" : "./",
                 "remote_profile":{
                     "email": "",
@@ -122,13 +122,13 @@
         },
         "bohrium5": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
             "machine_dict":{
-                "batch_type": "DpCloudServer",
-                "context_type": "DpCloudServerContext",
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
                 "local_root" : "./",
                 "remote_profile":{
                     "email": "",
@@ -161,6 +161,7 @@
         "run_select_config": "bohrium4",
         "prep_data_config": "bohrium3",
         "run_train_config": "bohrium5",
+        "model_devi_config": "bohrium5",
         "workflow_steps_config": "bohrium3"
     }
 }
\ No newline at end of file
diff --git a/rid/template/machine_bohrium_k8s.json b/rid/template/machine_bohrium_k8s.json
index b4378035..8fc03b57 100755
--- a/rid/template/machine_bohrium_k8s.json
+++ b/rid/template/machine_bohrium_k8s.json
@@ -4,18 +4,18 @@
             "template_config" : {
                 "image": "registry.dp.tech/public/pkufjhdocker/rid-tf-cpu:stable", 
                 "image_pull_policy" : "IfNotPresent",
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             }
         },
         "bohrium1": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
             "machine_dict":{
-                "batch_type": "DpCloudServer",
-                "context_type": "DpCloudServerContext",
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
                 "local_root" : "./",
                 "remote_profile":{
                     "email": "",
@@ -39,13 +39,13 @@
         },
         "bohrium2": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
             "machine_dict":{
-                "batch_type": "DpCloudServer",
-                "context_type": "DpCloudServerContext",
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
                 "local_root" : "./",
                 "remote_profile":{
                     "email": "",
@@ -69,13 +69,13 @@
         },
         "bohrium3": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
             "machine_dict":{
-                "batch_type": "DpCloudServer",
-                "context_type": "DpCloudServerContext",
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
                 "local_root" : "./",
                 "remote_profile":{
                     "email": "",
@@ -99,13 +99,13 @@
         },
         "bohrium4": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
             "machine_dict":{
-                "batch_type": "DpCloudServer",
-                "context_type": "DpCloudServerContext",
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
                 "local_root" : "./",
                 "remote_profile":{
                     "email": "",
@@ -138,6 +138,7 @@
         "run_select_config": "bohrium3",
         "prep_data_config": "local_k8s",
         "run_train_config": "bohrium4",
+        "model_devi_config": "bohrium4",
         "workflow_steps_config": "local_k8s"
     }
 }
\ No newline at end of file
diff --git a/rid/template/machine_bohrium_k8s_mcmc.json b/rid/template/machine_bohrium_k8s_mcmc.json
new file mode 100644
index 00000000..2c4b90a0
--- /dev/null
+++ b/rid/template/machine_bohrium_k8s_mcmc.json
@@ -0,0 +1,79 @@
+{
+    "resources": {
+        "local_k8s": {
+            "template_config" : {
+                "image": "registry.dp.tech/public/pkufjhdocker/rid-tf-cpu:stable", 
+                "image_pull_policy" : "IfNotPresent",
+                "requests" : {"ephemeral-storage": "4Gi"}
+            }
+        },
+        "bohrium1": {
+            "template_config":{
+                "requests" : {"ephemeral-storage": "4Gi"}
+            },
+            "executor":{
+            "image": "dptechnology/dpdispatcher:latest",
+            "machine_dict":{
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
+                "local_root" : "./",
+                "remote_profile":{
+                    "email": "",
+                    "password": "",
+                    "program_id": "",
+                    "input_data":{
+                        "api_version":2,
+                        "job_type": "container",
+                        "job_group_id": "",
+                        "log_file": "tmp_log",
+                        "job_name": "train",
+                        "scass_type":"c12_m92_1 * NVIDIA V100",
+                        "platform": "ali",
+                        "image_name":"registry.dp.tech/public/pkufjhdocker/rid-tf-gpu:stable"
+                        }
+                }
+            },
+            "resources_dict":{
+                "batch_type:": "Bohrium"
+            }
+        }
+        },
+        "bohrium2": {
+            "template_config":{
+                "requests" : {"ephemeral-storage": "4Gi"}
+            },
+            "executor":{
+            "image": "dptechnology/dpdispatcher:latest",
+            "machine_dict":{
+                "batch_type": "Bohrium",
+                "context_type": "Bohrium",
+                "local_root" : "./",
+                "remote_profile":{
+                    "email": "",
+                    "password": "",
+                    "program_id": "",
+                    "input_data":{
+                        "api_version":2,
+                        "job_type": "container",
+                        "job_group_id": "",
+                        "log_file": "tmp_log",
+                        "job_name": "train",
+                        "scass_type":"c12_m92_1 * NVIDIA V100",
+                        "platform": "ali",
+                        "image_name":"registry.dp.tech/public/pkufjhdocker/rid-tf-cpu:stable"
+                        }
+                }
+            },
+            "resources_dict":{
+                "batch_type:": "Bohrium"
+            }
+        }
+        }
+    },
+
+    "tasks": {
+        "mcmc_run_config": "bohrium1",
+        "mcmc_plot_config": "bohrium2",
+        "workflow_steps_config": "local_k8s"
+    }
+}
\ No newline at end of file
diff --git a/rid/template/machine_local_k8s.json b/rid/template/machine_local_k8s.json
index e8dbefe2..43b628d1 100644
--- a/rid/template/machine_local_k8s.json
+++ b/rid/template/machine_local_k8s.json
@@ -46,6 +46,7 @@
         "run_select_config": "local_k8s_4",
         "prep_data_config": "local_k8s_3",
         "run_train_config": "local_k8s_5",
+        "model_devi_config": "local_k8s_5",
         "workflow_steps_config": "local_k8s_3"
     }
 }
\ No newline at end of file
diff --git a/rid/template/machine_slurm_k8s.json b/rid/template/machine_slurm_k8s.json
index b23d36e7..b55efefa 100644
--- a/rid/template/machine_slurm_k8s.json
+++ b/rid/template/machine_slurm_k8s.json
@@ -4,12 +4,12 @@
             "template_config" : {
                 "image": "registry.dp.tech/public/pkufjhdocker/rid-tf-cpu:stable", 
                 "image_pull_policy" : "IfNotPresent",
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             }
         },
         "remote_slurm": {
             "template_config":{
-                "requests" : {"ephemeral-storage": "1Gi"}
+                "requests" : {"ephemeral-storage": "4Gi"}
             },
             "executor":{
             "image": "dptechnology/dpdispatcher:latest",
@@ -51,6 +51,7 @@
         "run_select_config": "remote_slurm",
         "prep_data_config": "local_k8s",
         "run_train_config": "remote_slurm",
+        "model_devi_config": "remote_slurm",
         "workflow_steps_config": "local_k8s"
     }
 }
\ No newline at end of file
diff --git a/rid/template/machine_slurm_k8s_mcmc.json b/rid/template/machine_slurm_k8s_mcmc.json
new file mode 100644
index 00000000..8be5069c
--- /dev/null
+++ b/rid/template/machine_slurm_k8s_mcmc.json
@@ -0,0 +1,50 @@
+{
+    "resources": {
+        "local_k8s": {
+            "template_config" : {
+                "image": "registry.dp.tech/public/pkufjhdocker/rid-tf-cpu:stable", 
+                "image_pull_policy" : "IfNotPresent",
+                "requests" : {"ephemeral-storage": "4Gi"}
+            }
+        },
+        "remote_slurm": {
+            "template_config":{
+                "requests" : {"ephemeral-storage": "4Gi"}
+            },
+            "executor":{
+            "image": "dptechnology/dpdispatcher:latest",
+            "machine_dict":{
+                "batch_type": "Slurm",
+                "context_type": "SSHContext",
+                "local_root" : "./",
+                "remote_root": "",
+                "remote_profile":{
+                    "hostname": "",
+                    "username": "",
+                    "password": "",
+                    "port": "",
+                    "timeout": 20
+                }
+            },
+            "resources_dict":{
+                "number_node": 1,
+                "cpu_per_node": 8,
+                "gpu_per_node": 1,
+                "queue_name": "",
+                "group_size": 1,
+                "custom_flags": [
+                    "#SBATCH --time=120:00:00",
+                    "#SBATCH --exclude=3090-[000,007]"
+                   ],
+                "source_list": ["/path/to/rid.env"]
+            }
+        }
+        }
+    },
+
+    "tasks": {
+        "mcmc_run_config": "remote_slurm",
+        "mcmc_plot_config": "remote_slurm",
+        "workflow_steps_config": "local_k8s"
+    }
+}
\ No newline at end of file
diff --git a/rid/template/machine_slurm_local.json b/rid/template/machine_slurm_local.json
index fc29c60e..9c22bd87 100644
--- a/rid/template/machine_slurm_local.json
+++ b/rid/template/machine_slurm_local.json
@@ -57,6 +57,7 @@
         "run_select_config": "local_machine",
         "prep_data_config": "local_machine",
         "run_train_config": "local_machine",
+        "model_devi_config": "local_machine",
         "workflow_steps_config": "local_machine_cpu"
     }
 }
\ No newline at end of file
diff --git a/rid/template/machine_slurm_local_mcmc.json b/rid/template/machine_slurm_local_mcmc.json
new file mode 100644
index 00000000..8e5c9f6a
--- /dev/null
+++ b/rid/template/machine_slurm_local_mcmc.json
@@ -0,0 +1,56 @@
+{
+    "resources": {
+        "local_machine": {
+            "executor":{
+            "image": "dptechnology/dpdispatcher:latest",
+            "merge_sliced_step": "True",
+            "machine_dict":{
+                "batch_type": "Slurm",
+                "context_type": "Local",
+                "local_root" : "./",
+                "remote_root": ""
+            },
+            "resources_dict":{
+                "number_node": 1,
+                "cpu_per_node": 8,
+                "gpu_per_node": 1,
+                "queue_name": "GPU",
+                "group_size": 1,
+                "custom_flags": [
+                    "#SBATCH --time=120:00:00"
+                   ],
+                "source_list": ["/path/to/rid.env"]
+            }
+        }
+        },
+        "local_machine_cpu": {
+            "executor":{
+            "image": "dptechnology/dpdispatcher:latest",
+            "merge_sliced_step": "True",
+            "machine_dict":{
+                "batch_type": "Slurm",
+                "context_type": "Local",
+                "local_root" : "./",
+                "remote_root": ""
+            },
+            "resources_dict":{
+                "number_node": 1,
+                "cpu_per_node": 8,
+                "gpu_per_node": 0,
+                "queue_name": "C032M0128G",
+                "group_size": 1,
+                "custom_flags": [
+                    "#SBATCH --time=120:00:00"
+                   ],
+                "source_list": ["/path/to/rid.env"]
+            }
+        }
+        }
+    },
+
+    "tasks": {
+        "mcmc_run_config": "remote_slurm",
+        "mcmc_plot_config": "remote_slurm",
+        "workflow_steps_config": "local_k8s"
+    }
+}
\ No newline at end of file
diff --git a/rid/template/rid_lmp_gmx_dp.json b/rid/template/rid_lmp_gmx_dp.json
new file mode 100755
index 00000000..d58c5a9c
--- /dev/null
+++ b/rid/template/rid_lmp_gmx_dp.json
@@ -0,0 +1,100 @@
+{
+    "name": "rid-mixed-gmx-lmp",
+    "numb_walkers": 2,
+    "numb_iters": 2,
+    "trust_lvl_1": 2,
+    "trust_lvl_2": 4,
+    "init_models": [],
+    
+    "CV": {
+        "mode": "custom",
+        "selected_atomid":[[161,165],[124,156],[124,161],[156,165],[107,161],[139,161],[153,161],[121,162],[121,163],[121,164],[138,161],[154,161],[153,165]],
+        "angular_mask": [0,0,0,0,0,0,0,0,0,0,0,0,0],
+        "weights": [1,1,1,1,1,1,1,1,1,1,1,1,1],
+        "units": "A",
+        "cv_file": ["colvar", "plmpath.pdb"]
+    },
+
+    "SelectorConfig": {
+        "numb_cluster_lower": 6,
+        "numb_cluster_upper": 10,
+        "cluster_threshold": 2.0,
+        "max_selection": 20,
+        "numb_cluster_threshold": 5,
+        "slice_mode": "dpdata",
+        "type_map": ["Si", "O", "H", "Hf", "Cl"]
+    },
+
+    "ExploreMDConfig": {
+        "type": "lmp",
+        "temperature": 372,
+        "dp_model": ["1a0.pb","1a1.pb","1a2.pb","1a3.pb"],
+        "dp_model_devi_out": "model_devi.out",
+        "dp_model_devi_f_lo": 0.14,
+        "dp_model_devi_f_hi": 0.26,
+        "dp_model_max_selection": 20,
+        "proj_info": {
+            "proj_mode": "path",
+            "proj_cv_index": [0,1,2,3],
+            "path_list": [[2.38722,0.99005,3.99886,2.08580],[2.64067,1.20917,3.73831,1.86004],[2.90639,1.50430,3.47935,1.65185],
+            [3.18130,1.83919,3.22240,1.46871],[3.46322,2.19574,2.96795,1.32108],[3.75057,2.56491,2.71672,1.22191],
+            [4.04218,2.94198,2.46968,1.18344],[4.33720,3.32424,2.22824,1.21147],[4.63498,3.71010,1.99443,1.30170]],
+            "path_lm": 0.25
+        },
+        "inputfile":"input_explore.lammps",
+        "dt": 0.001,
+        "nsteps": 10000,
+        "output_freq": 100,
+        "ntmpi": 1,
+        "nt": 12,
+        "max_warning": 2
+    },
+
+    "LabelMDConfig": {
+        "type": "gmx",
+        "method": "constrained",
+        "temperature": 372,
+        "dp_model": ["1a0.pb"],
+        "nsteps": 10000,
+        "ref-t": "372 0",
+        "gen-vel": "yes",
+        "gen-temp": "372",
+        "tc-grps": "top bottom",
+        "dt": 0.001,
+        "tau-t": "0.05 0.05",
+        "constraints": "none", 
+        "comm-grps": "top",
+        "output_freq": 100,
+        "output_mode": "both",
+        "freezegrps": "bottom",
+        "freezedim": "Y Y Y",
+        "pcoupl": "no",
+        "ntmpi": 1,
+        "coulombtype": "cutoff",
+        "rlist": 1,
+        "rvdw": 0.7,
+        "rvdw-switch": 0,
+        "rcoulomb": 0.7,
+        "rcoulomb-switch": 0,
+        "nstenergy": 1000,
+        "nt": 12,
+        "max_warning": 3
+    },
+
+    "Train": {
+        "numb_models": 4,
+        "neurons": [120, 120, 120, 120],
+        "resnet": true,
+        "batch_size": 4096,
+        "epoches": 30000,
+        "init_lr": 0.001,
+        "decay_steps": 200,
+        "decay_rate": 0.96,
+        "train_thread": 8,
+        "drop_out_rate": 0.1,
+        "numb_threads": 8,
+        "use_mix": false,
+        "restart": false,
+        "decay_steps_inner": 200
+    }
+}
\ No newline at end of file
diff --git a/rid/template/rid_mcmc_cv_dih.json b/rid/template/rid_mcmc_cv_dih.json
new file mode 100755
index 00000000..e2cbc6a1
--- /dev/null
+++ b/rid/template/rid_mcmc_cv_dih.json
@@ -0,0 +1,18 @@
+{
+        "name": "rid-mcmc",
+        "init_models": ["model_000.pb","model_001.pb","model_002.pb","model_003.pb"],
+        
+        "MCMC_Config": {
+            "cv_dimension": 2,
+            "numb_steps": 10000,
+            "numb_walkers": 500,
+            "temperature": 300,
+            "proj_info": {
+                "proj_mode": "cv",
+                "proj_cv_index": [0,1]
+            },
+            "cv_type": "dih",
+            "bins": 101
+        }
+    
+}
\ No newline at end of file
diff --git a/rid/template/rid_mcmc_cv_dis.json b/rid/template/rid_mcmc_cv_dis.json
new file mode 100755
index 00000000..5f0a707d
--- /dev/null
+++ b/rid/template/rid_mcmc_cv_dis.json
@@ -0,0 +1,20 @@
+{
+        "name": "rid-mcmc",
+        "init_models": ["model_000.pb","model_001.pb","model_002.pb","model_003.pb"],
+        
+        "MCMC_Config": {
+            "cv_dimension": 20,
+            "numb_steps": 100000,
+            "numb_walkers": 2000,
+            "temperature": 300,
+            "cv_upper_bound": [4.5,5.0,6.0,4.2,4.0,4.0,3.5,4.0,4.0,3.5,3.2,3.0,3.5,3.5,4.0,4.5,4.5,3.75,3.8,4.0],
+            "cv_lower_bound": [1.0,2.0,3.5,2.8,2.5,2.0,1.2,1.5,1.5,1.0,1.8,0.5,1.0,1.0,1.5,1.5,1.5,1.75,2.2,2.6],
+            "proj_info": {
+                "proj_mode": "cv",
+                "proj_cv_index": [0,1]
+            },
+            "cv_type": "dis",
+            "bins": 101
+        }
+    
+}
\ No newline at end of file
diff --git a/rid/template/rid_mcmc_path.json b/rid/template/rid_mcmc_path.json
new file mode 100755
index 00000000..451ea58f
--- /dev/null
+++ b/rid/template/rid_mcmc_path.json
@@ -0,0 +1,24 @@
+{
+        "name": "rid-mcmc",
+        "init_models": ["model_000.pb","model_001.pb","model_002.pb","model_003.pb"],
+        
+        "MCMC_Config": {
+            "cv_dimension": 13,
+            "numb_steps": 100000,
+            "numb_walkers": 2000,
+            "temperature": 300,
+            "cv_upper_bound": [5.33507,4.43402,4.59772,3.40287,4.67268,5.12813,3.28143,3.74967,3.58144,5.34143,7.13014,4.58024,7.41678],
+            "cv_lower_bound": [2.38722,0.99005,3.99886,2.08580,3.81799,3.93331,2.768175,3.145618,3.161886,4.661101,5.286805,4.076754,2.893047],
+            "proj_info": {
+                "proj_mode": "path",
+                "proj_cv_index": [0,1,2,3],
+                "path_list": [[2.38722,0.99005,3.99886,2.08580],[2.64067,1.20917,3.73831,1.86004],[2.90639,1.50430,3.47935,1.65185],
+                [3.18130,1.83919,3.22240,1.46871],[3.46322,2.19574,2.96795,1.32108],[3.75057,2.56491,2.71672,1.22191],
+                [4.04218,2.94198,2.46968,1.18344],[4.33720,3.32424,2.22824,1.21147],[4.63498,3.71010,1.99443,1.30170]],
+                "path_lm": 0.25
+            },
+            "cv_type": "dis",
+            "bins": 101
+        }
+    
+}
\ No newline at end of file
diff --git a/setup.py b/setup.py
index 580eaa3a..5acda4d9 100644
--- a/setup.py
+++ b/setup.py
@@ -36,6 +36,7 @@
         "rid/nn",
         "rid/op",
         "rid/select",
+        "rid/mcmc",
         "rid/superop",
         "rid/task",
         "rid/utils",


Test patch:
diff --git a/tests/data/conf.json b/tests/data/conf.json
new file mode 100644
index 00000000..03eeec12
--- /dev/null
+++ b/tests/data/conf.json
@@ -0,0 +1,1 @@
+{"conf.gro":"000"}
\ No newline at end of file
diff --git a/tests/data/model_000.pb b/tests/data/model_000.pb
deleted file mode 100644
index eedf9dc3..00000000
Binary files a/tests/data/model_000.pb and /dev/null differ
diff --git a/tests/data/model_001.pb b/tests/data/model_001.pb
deleted file mode 100644
index a06c0362..00000000
Binary files a/tests/data/model_001.pb and /dev/null differ
diff --git a/tests/data/model_002.pb b/tests/data/model_002.pb
deleted file mode 100644
index 368902ff..00000000
Binary files a/tests/data/model_002.pb and /dev/null differ
diff --git a/tests/data/models/model_000.pb b/tests/data/models/model_000.pb
new file mode 100644
index 00000000..1ed33b56
Binary files /dev/null and b/tests/data/models/model_000.pb differ
diff --git a/tests/data/models/model_001.pb b/tests/data/models/model_001.pb
new file mode 100644
index 00000000..f44b3b5d
Binary files /dev/null and b/tests/data/models/model_001.pb differ
diff --git a/tests/data/models/model_002.pb b/tests/data/models/model_002.pb
new file mode 100644
index 00000000..f439d618
Binary files /dev/null and b/tests/data/models/model_002.pb differ
diff --git a/tests/data/models/model_003.pb b/tests/data/models/model_003.pb
new file mode 100644
index 00000000..ca409ddf
Binary files /dev/null and b/tests/data/models/model_003.pb differ
diff --git a/tests/data/rid.json b/tests/data/rid.json
index 10279188..efafbc9f 100755
--- a/tests/data/rid.json
+++ b/tests/data/rid.json
@@ -25,6 +25,7 @@
 
     "ExploreMDConfig": {
         "nsteps": 50,
+        "type": "gmx",
         "output_freq": 1,
         "temperature": 300,
         "dt": 0.002,
@@ -36,6 +37,7 @@
 
     "LabelMDConfig": {
         "nsteps": 50,
+        "type": "gmx",
         "output_freq": 1,
         "temperature": 300,
         "dt": 0.002,
diff --git a/tests/mocked_ops.py b/tests/mocked_ops.py
index 686edb0c..131c7575 100755
--- a/tests/mocked_ops.py
+++ b/tests/mocked_ops.py
@@ -28,7 +28,7 @@
 from rid.op.run_train import TrainModel
 from rid.constants import (
     tf_model_name,
-    init_conf_name,
+    init_conf_gmx_name,
     data_new,
     data_raw,
     gmx_conf_name,
@@ -59,9 +59,13 @@ def clear_dir(data=None):
 
 def clear_files(files=None):
     if files:
-        for file in files:
-            if os.path.exists(file):
-                os.remove(file)
+        if isinstance(files,list):
+            for file in files:
+                if os.path.exists(file):
+                    os.remove(file)
+        elif isinstance(files,Path):
+            if os.path.exists(files):
+                    os.remove(files)
 
 def make_mocked_init_models(numb_models):
     tmp_models = []
@@ -85,7 +89,7 @@ def make_mocked_init_data(data_name,data_num):
 def make_mocked_init_confs(numb_confs):
     tmp_confs = []
     for ii in range(numb_confs):
-        ff = Path(init_conf_name.format(idx=ii))
+        ff = Path(init_conf_gmx_name.format(idx=ii))
         ff.write_text(f'This is init conf {ii}')
         tmp_confs.append(ff)
     return tmp_confs
@@ -95,6 +99,12 @@ def make_mocked_file(file_name):
     ff.write_text(f'This is mocked %s'%file_name)
     return ff
 
+def make_mocked_json(filename,fconts):
+    ff = Path(filename)
+    with open(ff,"w") as f:
+        json.dump(fconts,f)
+    return ff
+
 class MockedCollectData(CollectData):
     @OP.exec_sign_check
     def execute(
@@ -164,6 +174,7 @@ def execute(
             self,
             ip : OPIO,
     ) -> OPIO:
+        dp_conf_list = []
         assert(ip["task_path"].is_dir())
         print("start run explore !!")
         with set_directory(ip["task_path"]):
@@ -183,6 +194,7 @@ def execute(
                 f.write("gmx trajectory file")
         
         op = OPIO({
+                "dp_selected_confs": dp_conf_list,
                 "plm_out": ip["task_path"].joinpath(plumed_output_name),
                 "md_log": ip["task_path"].joinpath(gmx_mdrun_log),
                 "trajectory": ip["task_path"].joinpath(gmx_xtc_name),
@@ -204,12 +216,12 @@ def execute(
 
             tags = {}
             for tag in ip["conf_tags"]:
-                if isinstance(tag, Dict):
-                    tags.update(tag)
-                elif isinstance(tag, str):
-                    tags.update(json.loads(tag))
+                if isinstance(tag,Path):
+                    with open(tag,"r") as f:
+                        tags.update(json.load(f))
                 else:
                     raise RuntimeError("Unkown Error.")
+                
             conf_tags = []
             for conf in ip["confs"]:
                 conf_tags.append(str(tags[conf.name]))
@@ -323,6 +335,8 @@ def execute(
                 conf_tags[sel_gro_name.format(walker = int(ip['task_name']),idx=sel)] = f"{ip['task_name']}_{sel}"
                 save_txt(cv_init_label.format(walker = int(ip['task_name']),idx=sel), sel_data[ii])
                 cv_init_list.append(task_path.joinpath(cv_init_label.format(walker = int(ip['task_name']),idx=sel)))
+        with open("conf.json", "w") as f:
+            json.dump(conf_tags,f)
 
         op_out = OPIO(
             {
@@ -330,7 +344,7 @@ def execute(
                "selected_cv_init": cv_init_list,
                "model_devi": task_path.joinpath("cls_"+model_devi_name),
                "selected_indices": task_path.joinpath(sel_ndx_name),
-               "selected_conf_tags": conf_tags
+               "selected_conf_tags": task_path.joinpath("conf.json")
             }
         )
         return op_out
diff --git a/tests/op/test_prep_label.py b/tests/op/test_prep_label.py
index 8d770548..b344c1c0 100644
--- a/tests/op/test_prep_label.py
+++ b/tests/op/test_prep_label.py
@@ -14,14 +14,7 @@
 from pathlib import Path
 import shutil
 from rid.constants import (
-        explore_task_pattern, 
-        gmx_conf_name,
-        gmx_top_name,
-        gmx_mdp_name, 
-        plumed_input_name,
-        plumed_output_name,
-        sel_gro_name,
-        init_conf_name
+        plumed_input_name
     )
 class Test_CheckLabelInputs(unittest.TestCase):
     def setUp(self) -> None:
@@ -31,11 +24,12 @@ def test(self):
         op = CheckLabelInputs()
         data = Path(self.datapath)
         conf_path = data/"conf.gro"
+        conf_tags_path = data/"conf.json"
         
         op_in1 = OPIO(
             {
                 "confs": [conf_path],
-                "conf_tags": [{"conf.gro": "000"}]
+                "conf_tags": [conf_tags_path]
             }
         )
         op_in2 = OPIO(
diff --git a/tests/op/test_prep_rid.py b/tests/op/test_prep_rid.py
index 406b4f86..7e8015da 100644
--- a/tests/op/test_prep_rid.py
+++ b/tests/op/test_prep_rid.py
@@ -14,14 +14,7 @@
 from pathlib import Path
 import shutil
 from rid.constants import (
-        explore_task_pattern, 
-        gmx_conf_name,
-        gmx_top_name,
-        gmx_mdp_name, 
-        plumed_input_name,
-        plumed_output_name,
-        sel_gro_name,
-        init_conf_name
+        init_conf_gmx_name
     )
 
 class Test_PrepRid(unittest.TestCase):
@@ -34,7 +27,7 @@ def tearDown(self):
         data = Path(self.datapath)
         ii=data/self.confs
         for idx in range(self.numb_walkers):
-            os.remove(init_conf_name.format(idx=idx))
+            os.remove(init_conf_gmx_name.format(idx=idx))
     
     def test(self):
         op = PrepRiD()
diff --git a/tests/op/test_prep_select.py b/tests/op/test_prep_select.py
index 32dba677..e6e0b9cc 100644
--- a/tests/op/test_prep_select.py
+++ b/tests/op/test_prep_select.py
@@ -13,16 +13,6 @@
 from rid.utils import load_txt, save_txt, set_directory
 from pathlib import Path
 import shutil
-from rid.constants import (
-        explore_task_pattern, 
-        gmx_conf_name,
-        gmx_top_name,
-        gmx_mdp_name, 
-        plumed_input_name,
-        plumed_output_name,
-        sel_gro_name,
-        init_conf_name
-    )
 
 class Test_PrepSelect(unittest.TestCase):
     def setUp(self):
diff --git a/tests/op/test_run_label.py b/tests/op/test_run_label.py
index b1bd3c21..84ec3f4d 100644
--- a/tests/op/test_run_label.py
+++ b/tests/op/test_run_label.py
@@ -14,16 +14,6 @@
 from rid.utils import load_txt, save_txt, set_directory
 from pathlib import Path
 import shutil
-from rid.constants import (
-        explore_task_pattern, 
-        gmx_conf_name,
-        gmx_top_name,
-        gmx_mdp_name, 
-        plumed_input_name,
-        plumed_output_name,
-        sel_gro_name,
-        init_conf_name
-    )
 
 class Test_MockedRunLabel(unittest.TestCase):
     def setUp(self):
diff --git a/tests/op/test_run_select.py b/tests/op/test_run_select.py
index 975ee584..b6a2b4a9 100644
--- a/tests/op/test_run_select.py
+++ b/tests/op/test_run_select.py
@@ -14,16 +14,6 @@
 from rid.utils import load_txt, save_txt, set_directory
 from pathlib import Path
 import shutil
-from rid.constants import (
-        explore_task_pattern, 
-        gmx_conf_name,
-        gmx_top_name,
-        gmx_mdp_name, 
-        plumed_input_name,
-        plumed_output_name,
-        sel_gro_name,
-        init_conf_name
-    )
 
 class Test_MockedRunSelect(unittest.TestCase):
     def setUp(self):
@@ -53,9 +43,11 @@ def test(self, mocked_run):
         cluster_data = data/"cls_sel.out.npy"
         xtc_traj = data/"traj_comp.xtc"
         top_path = data/"topol.top"
-        model_path1 = ".."/data/"model_000.pb"
-        model_path2 = ".."/data/"model_001.pb"
-        model_path3 = ".."/data/"model_002.pb"
+        model_path1 = ".."/data/"models/model_000.pb"
+        model_path2 = ".."/data/"models/model_001.pb"
+        model_path3 = ".."/data/"models/model_002.pb"
+        gmx_config = {"type":"gmx","nsteps": 50,"method":"restrained", "output_freq": 1, "temperature": 300, "kappas": [500,500],
+                      "dt": 0.002, "output_mode": "both", "ntmpi": 1, "nt": 8, "max_warning": 0}
 
         op_in1 = OPIO(
             {
@@ -67,6 +59,7 @@ def test(self, mocked_run):
                 "trust_lvl_2": 0.03,
                 "xtc_traj": xtc_traj,
                 "topology": top_path,
+                "label_config": gmx_config,
                 "dt": 0.002,
                 "output_freq": 2500,
                 "slice_mode": "gmx"
@@ -82,6 +75,7 @@ def test(self, mocked_run):
                 "trust_lvl_2": 0.03,
                 "xtc_traj": xtc_traj,
                 "topology": top_path,
+                "label_config": gmx_config,
                 "dt": 0.002,
                 "output_freq": 2500,
                 "slice_mode": "mdtraj"
@@ -97,6 +91,7 @@ def test(self, mocked_run):
                 "trust_lvl_2": 0.03,
                 "xtc_traj": xtc_traj,
                 "topology": top_path,
+                "label_config": gmx_config,
                 "dt": 0.002,
                 "output_freq": 2500,
                 "slice_mode": "gmx"
@@ -112,6 +107,7 @@ def test(self, mocked_run):
                 "trust_lvl_2": 0.03,
                 "xtc_traj": xtc_traj,
                 "topology": top_path,
+                "label_config": gmx_config,
                 "dt": 0.002,
                 "output_freq": 2500,
                 "slice_mode": "others"
diff --git a/tests/test_label.py b/tests/test_label.py
index 43aef1c1..9122b83d 100755
--- a/tests/test_label.py
+++ b/tests/test_label.py
@@ -39,6 +39,7 @@
     make_mocked_init_data,
     make_mocked_init_models,
     make_mocked_init_confs,
+    make_mocked_json,
     MockedCheckLabelInputs,
     MockedPrepLabel,
     MockedRunLabel
@@ -140,7 +141,8 @@ def setUp(self):
         self.kappas = [500.]*3
         self.angular_mask = [0,0,0]
         self.tail = 0.9
-        self.conf_tags = [{"conf_000.gro":"000_000"},{"conf_001.gro":"001_001"},{"conf_002.gro":"002_002"}]
+        self.init_conf_tags = make_mocked_json("conf.json", {"conf_000.gro":"000_000", "conf_001.gro":"001_001", "conf_002.gro":"002_002"})
+        self.conf_tags = upload_artifact(self.init_conf_tags)
         self.label_config = label_config
         self.cv_config = cv_config
         self.data_out = "label"
@@ -149,6 +151,7 @@ def setUp(self):
 
     def tearDown(self):
         clear_files(self.init_models)
+        clear_files(self.init_conf_tags)
         clear_files(self.init_confs)
         clear_dir(self.center[0])
         clear_dir(self.data_out)
@@ -169,10 +172,10 @@ def test_label(self):
                 "label_config" : self.label_config,
                 "cv_config" : self.cv_config,
                 "tail": self.tail,
-                "conf_tags": self.conf_tags,
                 "block_tag" : self.block_tag
             },
             artifacts = {
+                "conf_tags": self.conf_tags,
                 "models" : None,
                 "forcefield": None,
                 "topology" : None,
diff --git a/tests/test_select.py b/tests/test_select.py
index 8c778ab6..d5f2071b 100644
--- a/tests/test_select.py
+++ b/tests/test_select.py
@@ -74,6 +74,7 @@ def setUp(self):
         self.task_name = ["000","001","002"]
         self.block_tag = "000"
         self.data_out = "select"
+        self.label_config = {"type":"gmx"}
 
     def tearDown(self):
         clear_dir(self.plm_out[0])
@@ -93,6 +94,7 @@ def test_select(self):
             'select-step',
             template = steps,
             parameters = {
+                "label_config": self.label_config,
                 "trust_lvl_1" : self.trust_lvl_1,
                 "trust_lvl_2": self.trust_lvl_2,
                 "cluster_threshold": self.cluster_threshold,