fixed swap

Author: simongravelle

Diff for: docs/source/chapters/chapter3.rst

@@ -183,8 +183,8 @@ neighboring atoms, the simulation becomes more efficient. Add the following
 
 .. code-block:: python
 
-    def update_neighbor_lists(self):
-        if (self.step % self.neighbor == 0):
+    def update_neighbor_lists(self, force_update=False):
+        if (self.step % self.neighbor == 0) | force_update:
             matrix = distances.contact_matrix(self.atoms_positions,
                 cutoff=self.cut_off, #+2,
                 returntype="numpy",
@@ -218,8 +218,8 @@ below).
 
 .. code-block:: python
 
-    def update_cross_coefficients(self):
-        if (self.step % self.neighbor == 0):
+    def update_cross_coefficients(self, force_update=False):
+        if (self.step % self.neighbor == 0) | force_update:
             # Precalculte LJ cross-coefficients
             sigma_ij_list = []
             epsilon_ij_list = []

Diff for: docs/source/chapters/chapter8.rst

@@ -29,39 +29,40 @@ Let us add the following method called *monte_carlo_exchange* to the *MonteCarlo
 .. code-block:: python
 
     def monte_carlo_exchange(self):
-        # The first step is to make a copy of the previous state
-        # Since atoms numbers are evolving, its also important to store the
-        # neighbor, sigma, and epsilon lists
-        self.Epot = self.compute_potential() # TOFIX: not necessary every time
-        initial_positions = copy.deepcopy(self.atoms_positions)
-        initial_number_atoms = copy.deepcopy(self.number_atoms)
-        initial_neighbor_lists = copy.deepcopy(self.neighbor_lists)
-        initial_sigma_lists = copy.deepcopy(self.sigma_ij_list)
-        initial_epsilon_lists = copy.deepcopy(self.epsilon_ij_list)
-        # Apply a 50-50 probability to insert or delete
-        insert_or_delete = np.random.random()
-        if np.random.random() < insert_or_delete:
-            self.monte_carlo_insert()
-        else:
-            self.monte_carlo_delete()
-        if np.random.random() < self.acceptation_probability: # accepted move
-            # Update the success counters
+        if self.desired_mu is not None:
+            # The first step is to make a copy of the previous state
+            # Since atoms numbers are evolving, its also important to store the
+            # neighbor, sigma, and epsilon lists
+            self.Epot = self.compute_potential() # TOFIX: not necessary every time
+            initial_positions = copy.deepcopy(self.atoms_positions)
+            initial_number_atoms = copy.deepcopy(self.number_atoms)
+            initial_neighbor_lists = copy.deepcopy(self.neighbor_lists)
+            initial_sigma_lists = copy.deepcopy(self.sigma_ij_list)
+            initial_epsilon_lists = copy.deepcopy(self.epsilon_ij_list)
+            # Apply a 50-50 probability to insert or delete
+            insert_or_delete = np.random.random()
             if np.random.random() < insert_or_delete:
-                self.successful_insert += 1
+                self.monte_carlo_insert()
             else:
-                self.successful_delete += 1
-        else:
-            # Reject the new position, revert to inital position
-            self.neighbor_lists = initial_neighbor_lists
-            self.sigma_ij_list = initial_sigma_lists
-            self.epsilon_ij_list = initial_epsilon_lists
-            self.atoms_positions = initial_positions
-            self.number_atoms = initial_number_atoms
-            # Update the failed counters
-            if np.random.random() < insert_or_delete:
-                self.failed_insert += 1
+                self.monte_carlo_delete()
+            if np.random.random() < self.acceptation_probability: # accepted move
+                # Update the success counters
+                if np.random.random() < insert_or_delete:
+                    self.successful_insert += 1
+                else:
+                    self.successful_delete += 1
             else:
-                self.failed_delete += 1
+                # Reject the new position, revert to inital position
+                self.neighbor_lists = initial_neighbor_lists
+                self.sigma_ij_list = initial_sigma_lists
+                self.epsilon_ij_list = initial_epsilon_lists
+                self.atoms_positions = initial_positions
+                self.number_atoms = initial_number_atoms
+                # Update the failed counters
+                if np.random.random() < insert_or_delete:
+                    self.failed_insert += 1
+                else:
+                    self.failed_delete += 1
 
 .. label:: end_MonteCarlo_class
 

Diff for: docs/source/chapters/chapter9.rst

@@ -36,6 +36,17 @@ Carlo swap move can considerably speed up equilibration.
             position2 = copy.deepcopy(self.atoms_positions[shift_2+atom_id_2])
             self.atoms_positions[shift_2+atom_id_2] = position1
             self.atoms_positions[shift_1+atom_id_1] = position2
+            # force the recalculation of neighbor list
+            initial_atoms_sigma = self.atoms_sigma
+            initial_atoms_epsilon = self.atoms_epsilon
+            initial_atoms_mass = self.atoms_mass
+            initial_atoms_type = self.atoms_type
+            initial_sigma_ij_list = self.sigma_ij_list
+            initial_epsilon_ij_list = self.epsilon_ij_list
+            initial_neighbor_lists = self.neighbor_lists
+            self.update_neighbor_lists(force_update=True)
+            self.identify_atom_properties()
+            self.update_cross_coefficients(force_update=True)
             # Measure the potential energy of the new configuration
             trial_Epot = self.compute_potential()
             # Evaluate whether the new configuration should be kept or not
@@ -49,7 +60,14 @@ Carlo swap move can considerably speed up equilibration.
             else: # Reject new position
                 self.atoms_positions = initial_positions # Revert to initial positions
                 self.failed_swap += 1
-
+                self.atoms_sigma = initial_atoms_sigma
+                self.atoms_epsilon = initial_atoms_epsilon
+                self.atoms_mass = initial_atoms_mass
+                self.atoms_type = initial_atoms_type
+                self.sigma_ij_list = initial_sigma_ij_list
+                self.epsilon_ij_list = initial_epsilon_ij_list
+                self.neighbor_lists = initial_neighbor_lists
+                
 .. label:: end_MonteCarlo_class
 
 Let us initialise swap counter: