Removed the YAML parameter acceptance and updated the documentation a…

…nd unit tests accordingly

docs/simulations.rst

@@ -279,9 +279,6 @@ include parameters for data analysis here.
   - :code:`proposal`: (Optional, Default: :code:`exhaustive`)
       The method for proposing simulations to be swapped. Available options include :code:`single`, :code:`neighboring`, and :code:`exhaustive`.
       For more details, please refer to :ref:`doc_proposal`.
-  - :code:`acceptance`: (Optional, Default: :code:`metropolis`)
-      The Monte Carlo method for swapping simulations. Available options include :code:`same-state`/:code:`same_state`, :code:`metropolis`, and :code:`metropolis-eq`/:code:`metropolis_eq`. 
-      For more details, please refer to :ref:`doc_acceptance`.
   - :code:`w_combine`: (Optional, Default: :code:`False`)
       Whether to perform weight combination or not. Note that weights averaged over from the last time the Wang-Landau incrementor was updated (instead of 
       final weights) will be used for weight combination. For more details about weight combination, please refer to :ref:`doc_w_schemes`.
@@ -375,7 +372,6 @@ infinity internally.
     nst_sim: null
     add_swappables: null
     proposal: 'exhaustive'
-    acceptance: 'metropolis' 
     w_combine: False
     w_mean_type: 'simple'
     N_cutoff: 1000

ensemble_md/replica_exchange_EE.py

@@ -116,13 +116,8 @@ def set_params(self, analysis):
 
               - If a required parameter is not specified in the YAML file.
               - If a specified parameter is not recognizable.
-              - If a specified weight combining scheme is not available.
-              - If a specified acceptance scheme is not available.
-              - If a specified free energy estimator is not available.
-              - If a specified method for error estimation is not available.
-              - If an integer parameter is not an integer.
-              - If a positive parameter is not positive.
-              - If a non-negative parameter is negative.
+              - If a specified option is not available for a parameter.
+              - If the data type or range (e.g., positive or negative) of a parameter is not correct.
               - If an invalid MDP file is detected.
         """
         self.warnings = []  # Store warnings, if any.
@@ -161,7 +156,6 @@ def set_params(self, analysis):
             "modify_coords": None,
             "nst_sim": None,
             "proposal": 'exhaustive',
-            "acceptance": "metropolis",
             "w_combine": False,
             "w_mean_type": 'simple',
             "N_cutoff": 1000,
@@ -197,9 +191,6 @@ def set_params(self, analysis):
         if self.proposal not in [None, 'single', 'neighboring', 'exhaustive']:  # deprecated option: multiple
             raise ParameterError("The specified proposal scheme is not available. Available options include 'single', 'neighboring', and 'exhaustive'.")  # noqa: E501
 
-        if self.acceptance not in [None, 'same-state', 'same_state', 'metropolis']:
-            raise ParameterError("The specified acceptance scheme is not available. Available options include 'same-state' and 'metropolis'.")  # noqa: E501
-
         if self.df_method not in [None, 'TI', 'BAR', 'MBAR']:
             raise ParameterError("The specified free energy estimator is not available. Available options include 'TI', 'BAR', and 'MBAR'.")  # noqa: E501
 
@@ -491,7 +482,6 @@ def print_params(self, params_analysis=False):
         print(f'Simulation inputs: {gro_str}, {top_str}, {self.mdp}')
         print(f"Verbose log file: {self.verbose}")
         print(f"Proposal scheme: {self.proposal}")
-        print(f"Acceptance scheme for swapping simulations: {self.acceptance}")
         print(f"Whether to perform weight combination: {self.w_combine}")
         print(f"Type of means for weight combination: {self.w_mean_type}")
         print(f"Whether to perform histogram correction: {self.hist_corr}")
@@ -999,45 +989,37 @@ def calc_prob_acc(self, swap, dhdl_files, states, shifts):
         prob_acc : float
             The acceptance ratio.
         """
-        if self.acceptance == "same_state" or self.acceptance == "same-state":
-            if states[swap[0]] == states[swap[1]]:  # same state, swap!
-                prob_acc = 1  # This must lead to an output of swap_bool = True from the function accept_or_reject
-            else:
-                prob_acc = 0  # This must lead to an output of swap_bool = False from the function accept_or_reject
-
-        else:  # i.e. must be metropolis, which requires the calculation of dU
-            # Now we calculate dU
-            if self.verbose is True:
-                print("  Proposing a move from (x^i_m, x^j_n) to (x^i_n, x^j_m) ...")
-            f0, f1 = dhdl_files[swap[0]], dhdl_files[swap[1]]
-            h0, h1 = get_headers(f0), get_headers(f1)
-            data_0, data_1 = (
-                extract_dataframe(f0, headers=h0).iloc[-1],
-                extract_dataframe(f1, headers=h1).iloc[-1],
+        if self.verbose is True:
+            print("  Proposing a move from (x^i_m, x^j_n) to (x^i_n, x^j_m) ...")
+        f0, f1 = dhdl_files[swap[0]], dhdl_files[swap[1]]
+        h0, h1 = get_headers(f0), get_headers(f1)
+        data_0, data_1 = (
+            extract_dataframe(f0, headers=h0).iloc[-1],
+            extract_dataframe(f1, headers=h1).iloc[-1],
+        )
+
+        # \Delta H to all states at the last time frame
+        # Notably, the can be regarded as H for each state since the reference state will have a value of 0 anyway.
+        dhdl_0 = data_0[-self.n_sub:]
+        dhdl_1 = data_1[-self.n_sub:]
+
+        # Old local index, which will only be used in "metropolis"
+        old_state_0 = states[swap[0]] - shifts[swap[0]]
+        old_state_1 = states[swap[1]] - shifts[swap[1]]
+
+        # New local index. Note that states are global indices, so we shift them back to the local ones
+        new_state_0 = states[swap[1]] - shifts[swap[0]]  # new state index (local index in simulation swap[0])
+        new_state_1 = states[swap[0]] - shifts[swap[1]]  # new state index (local index in simulation swap[1])
+
+        dU_0 = (dhdl_0[new_state_0] - dhdl_0[old_state_0]) / self.kT  # U^{i}_{n} - U^{i}_{m}, i.e. \Delta U (kT) to the new state  # noqa: E501
+        dU_1 = (dhdl_1[new_state_1] - dhdl_1[old_state_1]) / self.kT  # U^{j}_{m} - U^{j}_{n}, i.e. \Delta U (kT) to the new state  # noqa: E501
+        dU = dU_0 + dU_1
+        if self.verbose is True:
+            print(
+                f"  U^i_n - U^i_m = {dU_0:.2f} kT, U^j_m - U^j_n = {dU_1:.2f} kT, Total dU: {dU:.2f} kT"
             )
 
-            # \Delta H to all states at the last time frame
-            # Notably, the can be regarded as H for each state since the reference state will have a value of 0 anyway.
-            dhdl_0 = data_0[-self.n_sub:]
-            dhdl_1 = data_1[-self.n_sub:]
-
-            # Old local index, which will only be used in "metropolis"
-            old_state_0 = states[swap[0]] - shifts[swap[0]]
-            old_state_1 = states[swap[1]] - shifts[swap[1]]
-
-            # New local index. Note that states are global indices, so we shift them back to the local ones
-            new_state_0 = states[swap[1]] - shifts[swap[0]]  # new state index (local index in simulation swap[0])
-            new_state_1 = states[swap[0]] - shifts[swap[1]]  # new state index (local index in simulation swap[1])
-
-            dU_0 = (dhdl_0[new_state_0] - dhdl_0[old_state_0]) / self.kT  # U^{i}_{n} - U^{i}_{m}, i.e. \Delta U (kT) to the new state  # noqa: E501
-            dU_1 = (dhdl_1[new_state_1] - dhdl_1[old_state_1]) / self.kT  # U^{j}_{m} - U^{j}_{n}, i.e. \Delta U (kT) to the new state  # noqa: E501
-            dU = dU_0 + dU_1
-            if self.verbose is True:
-                print(
-                    f"  U^i_n - U^i_m = {dU_0:.2f} kT, U^j_m - U^j_n = {dU_1:.2f} kT, Total dU: {dU:.2f} kT"
-                )
-
-            prob_acc = min(1, np.exp(-dU))
+        prob_acc = min(1, np.exp(-dU))
 
         return prob_acc
 

ensemble_md/tests/test_replica_exchange_EE.py

@@ -83,7 +83,6 @@ def test_set_params_error(self, params_dict):
 
         # 2. Available options
         check_param_error(params_dict, 'proposal', "The specified proposal scheme is not available. Available options include 'single', 'neighboring', and 'exhaustive'.", 'cool', 'exhaustive')  # noqa: E501
-        check_param_error(params_dict, 'acceptance', "The specified acceptance scheme is not available. Available options include 'same-state' and 'metropolis'.")  # noqa: E501
         check_param_error(params_dict, 'df_method', "The specified free energy estimator is not available. Available options include 'TI', 'BAR', and 'MBAR'.")  # noqa: E501
         check_param_error(params_dict, 'err_method', "The specified method for error estimation is not available. Available options include 'propagate', and 'bootstrap'.")  # noqa: E501
 
@@ -166,7 +165,6 @@ def test_set_params(self, params_dict):
 
         # 2. Check the default values of the parameters not specified in params.yaml
         assert REXEE.proposal == "exhaustive"
-        assert REXEE.acceptance == "metropolis"
         assert REXEE.w_combine is False
         assert REXEE.N_cutoff == 1000
         assert REXEE.verbose is True
@@ -266,7 +264,6 @@ def test_print_params(self, capfd, params_dict):
         L += "Simulation inputs: ensemble_md/tests/data/sys.gro, ensemble_md/tests/data/sys.top, ensemble_md/tests/data/expanded.mdp\n"  # noqa: E501
         L += "Verbose log file: True\n"
         L += "Proposal scheme: exhaustive\n"
-        L += "Acceptance scheme for swapping simulations: metropolis\n"
         L += "Whether to perform weight combination: False\n"
         L += "Type of means for weight combination: simple\n"
         L += "Whether to perform histogram correction: False\n"
@@ -492,32 +489,20 @@ def test_calc_prob_acc(self, capfd, params_dict):
         shifts = [0, 1, 2, 3]
         dhdl_files = [os.path.join(input_path, f"dhdl/dhdl_{i}.xvg") for i in range(4)]
 
-        # Test 1: Same-state swapping (True)
-        swap = (1, 2)
-        REXEE.acceptance = "same_state"
-        prob_acc_1 = REXEE.calc_prob_acc(swap, dhdl_files, states, shifts)
-        assert prob_acc_1 == 1
-
-        # Test 2: Same-state swapping (False)
-        swap = (0, 2)
-        prob_acc_2 = REXEE.calc_prob_acc(swap, dhdl_files, states, shifts)
-        assert prob_acc_2 == 0
-
-        # Test 3-1: Metropolis, test 1
+        # Test 1
         swap = (0, 1)
-        REXEE.acceptance = "metropolis"
-        prob_acc_3_1 = REXEE.calc_prob_acc(swap, dhdl_files, states, shifts)
+        prob_acc_1 = REXEE.calc_prob_acc(swap, dhdl_files, states, shifts)
         out, err = capfd.readouterr()
         # dU = (-9.1366697  + 11.0623788)/2.4777098766670016 ~ 0.7772 kT, so p_acc = 0.45968522728859024
-        assert prob_acc_3_1 == pytest.approx(0.45968522728859024)
+        assert prob_acc_1 == pytest.approx(0.45968522728859024)
         assert 'U^i_n - U^i_m = -3.69 kT, U^j_m - U^j_n = 4.46 kT, Total dU: 0.78 kT' in out
 
-        # Test 3-2: Metropolis, test 2
+        # Test 2
         swap = (0, 2)
-        prob_acc_3_2 = REXEE.calc_prob_acc(swap, dhdl_files, states, shifts)
+        prob_acc_2 = REXEE.calc_prob_acc(swap, dhdl_files, states, shifts)
         out, err = capfd.readouterr()
         # dU = (-9.1366697 + 4.9963939)/2.4777098766670016 ~ -1.6710 kT, so p_acc = 1
-        assert prob_acc_3_2 == 1
+        assert prob_acc_2 == 1
 
     def test_accept_or_reject(self, params_dict):
         REXEE = get_REXEE_instance(params_dict)