diff --git a/onix/cell.py b/onix/cell.py
index d52a0f1..9019063 100644
--- a/onix/cell.py
+++ b/onix/cell.py
@@ -543,9 +543,14 @@ def _change_total_density(self, s):
 
 			if self.name in density_change_dict:
 				bucell_density_change_dict = density_change_dict[self.name]
-				if s in bucell_density_change_dict:
-
-					new_tot_dens = bucell_density_change_dict[s]
+				# When user sets a new density for step s+1, the Monte Carlo simulation of step s+1
+				# needs to run with this new density
+				# For conveniency, ONIX thus changes the density at the end of step s
+				# This means that if we are now at step s, we need to look if the user has specified a new 
+				# density for step s+1
+				if s+1 in bucell_density_change_dict:
+
+					new_tot_dens = bucell_density_change_dict[s+1]
 					current_total_dens = self.get_total_dens()
 					factor = new_tot_dens/current_total_dens
 					passport_list = self.passlist.passport_list
@@ -560,9 +565,19 @@ def _change_isotope_density(self,s):
 		if isotopic_change_dict != None:
 			if self.name in isotopic_change_dict:
 				bucell_isotopic_change_dict = isotopic_change_dict[self.name]
+				unit = bucell_isotopic_change_dict['unit']
+				current_total_dens = self.get_total_dens()
 				for nucl_name in bucell_isotopic_change_dict:
+					if nucl_name == 'unit':
+						continue
 					nuclide_isotopic_change_dict = bucell_isotopic_change_dict[nucl_name]
-					if s in nuclide_isotopic_change_dict:
+
+					# When user sets a new density for step s+1, the Monte Carlo simulation of step s+1
+					# needs to run with this new density
+					# For conveniency, ONIX thus changes the density at the end of step s
+					# This means that if we are now at step s, we need to look if the user has specified a new 
+					# density for step s+1
+					if s+1 in nuclide_isotopic_change_dict:
 						nucl_passport = self.get_nuclide(nucl_name)
 						# Only the current_dens is set to the user-defined density
 						# the dens_seq and dens_subseq_mat last value are left unchanged
@@ -570,7 +585,14 @@ def _change_isotope_density(self,s):
 						# and mat_builder (which prepares the matrix for new depletion calculation) uses
 						# current_dens. Therefore, only changing current_dens will update calculation without
 						# changing stored value that will be printed
-						new_dens = nuclide_isotopic_change_dict[s]
+
+						# If unit is number density
+						if unit == 'number density':
+							new_dens = nuclide_isotopic_change_dict[s+1]
+						# else if unit is atom fraction
+						if unit == 'atom fraction':
+							new_dens = nuclide_isotopic_change_dict[s+1]*current_total_dens
+
 						nucl_passport.current_dens = new_dens
 
 
diff --git a/onix/sequence.py b/onix/sequence.py
index d06751c..f661985 100644
--- a/onix/sequence.py
+++ b/onix/sequence.py
@@ -665,17 +665,26 @@ def norma_unit(self, norma_unit):
     def isotopic_change_dict(self):
         return self._isotopic_change_dict
 
-    def set_isotopic_change(self, cell, cell_isotopic_change):
-        """Manually changes the isotopic densities (in atm per cm3) of user-specified nuclides in the BUCell for user-specified
+    def set_isotopic_change(self, cell, cell_isotopic_change, unit='number density'):
+        """Manually changes the isotopic densities of user-specified nuclides in the BUCell for user-specified
         macrosteps
+
+        Parameter:
+
+        unit: specifies the unit of the new isotopic density
+        'number density' (default) for density in atm per cm3
+        'atom fraction' for density as atomic fraction
         
-        IMPORTANT: The other method "set_density_change" trumps this method, i.e., the new isotopic densities set in
-        set_isotopic_change will be renormalized by the new value from set_density_change"""
+        IMPORTANT: If the method "set_density_change" is also used for the same BUCell, the isotopic change set by the user
+        must be in atom fraction"""
         # If this is the first cell which is set isotopic change, create the dict
         if self.isotopic_change_dict == None:
             self._isotopic_change_dict = {}
 
         self._isotopic_change_dict[cell.name] = cell_isotopic_change
+        self._isotopic_change_dict[cell.name]['unit'] = unit
+
+
 
 
     @property