Merge branch 'anneal' into 'master'

Anneal

See merge request sajuhak/oxidation_sim!7

config/default.in

@@ -4,9 +4,12 @@ report_file = True
 
 [Simulation parameters]
 sim_type = ase_sim
-kT = 300
-n_exchanges = 10000
-n_frames = 100
+kT_high = 500
+kT_low = 100
+n_steps_high = 10000
+n_steps_fall = 10000
+n_steps_low = 10000
+n_frames = 200
 rand_ox = True
 
 [Lattice parameters]
@@ -31,3 +34,4 @@ Cr = o,#8a99c7
 Al = ^,#bfa6a6
 O = D,#ff0d0d
 X = x,#000000
+

ovito_overlay.py

+import ovito
+
+# This user-defined function is called by OVITO to let it draw arbitrary graphics on top of the viewport.
+# It is passed a QPainter (see http://qt-project.org/doc/qt-5/qpainter.html).
+def render(painter, **args):
+
+	# This demo code prints the current animation frame into the upper left corner of the viewport.
+	text1 = "Frame {}".format(ovito.dataset.anim.current_frame)
+	y = 10 + painter.fontMetrics().ascent()
+	painter.drawText(10, y, text1)
+
+	node = ovito.dataset.selected_node
+	y += 10 + painter.fontMetrics().ascent()
+	text2 = "Temperature: {}".format(node.source.attributes['T'])
+	painter.drawText(10, y, text2)
+
+	y += 10 + painter.fontMetrics().ascent()
+	text3 = "Energy: {}".format(node.source.attributes['energy'])
+	painter.drawText(10, y, text3)
+
+	# Also print the current number of particles into the lower left corner of the viewport.
+	num_particles = (node.compute().number_of_particles if node else 0)
+	text4 = "{} particles".format(num_particles)
+	painter.drawText(10, painter.window().height() - 10, text4)
+
+	# Print to the log window:
+	print(text1)
+	print(text2)
+	print(text3)
+	print(text4)

sajuhak_dyn.py

@@ -101,6 +101,9 @@ class SajuhakDynamics(Dynamics):
 
     def run(self, kT, n_steps=100):
         for i in range(n_steps):
+            if type(kT) is np.ndarray:
+                self.step(kT[i])
+            else:
                 self.step(kT)
 
     def step(self, kT):

sim_controller.py

@@ -26,8 +26,8 @@ class SimController(object):
                      "{lattice_counts}"
                      "Particles created in {particles_time}.\n"
                      "\n"
-                     "Metropolis steps: {n_exchanges}\n"
-                     "kT: {kT}\n"
+                     "{kT_str}\n"
+                     "{exchanges_str}\n"
                      "Frames: {n_frames}\n"
                      "Random (de)oxidation: {rand_ox}\n"
                      "Simulation took {sim_time}.\n")
@@ -58,11 +58,30 @@ class SimController(object):
             sim_params = default_reader['Simulation parameters']
         if sim_params['sim_type'] == 'ase_sim':
             self.sim = ase_sim.AseSim()
-        self.n_exchanges = int(sim_params['n_exchanges'])
         self.n_frames = int(sim_params['n_frames'])
-        self.n_steps = int(self.n_exchanges / self.n_frames)
+        self.n_exchanges = None
+        self.kT = None
+        if 'n_exchanges' in sim_params.keys() and 'kT' in sim_params.keys():
+            self.n_exchanges = int(sim_params['n_exchanges'])
+            self.n_steps_per_frame = int(self.n_exchanges / self.n_frames)
             self.kT = int(sim_params['kT'])
-        if 'rand_ox' in sim_params.keys() and sim_params['rand_ox'] != False:
+            self.last_kT = self.kT
+        elif 'n_steps_high' in sim_params.keys() \
+                and 'n_steps_fall' in sim_params.keys() \
+                and 'n_steps_low' in sim_params.keys()\
+                and 'kT_high' in sim_params.keys()\
+                and 'kT_low' in sim_params.keys():
+            self.n_exchanges = [int(sim_params['n_steps_high']),
+                                 int(sim_params['n_steps_fall']),
+                                 int(sim_params['n_steps_low'])]
+            self.n_steps_per_frame = int(sum(self.n_exchanges) / self.n_frames)
+            self.kT = [int(sim_params['kT_high']),
+                                 int(sim_params['kT_low'])]
+            self.last_kT = self.kT[0]
+        else:
+            print('Could not determine the amount of steps.')
+            exit(1)
+        if 'rand_ox' in sim_params.keys() and sim_params['rand_ox'] is not False:
             self.rand_ox = True
         else:
             self.rand_ox = False
@@ -102,6 +121,8 @@ class SimController(object):
         else:
             self.figure_params = None
 
+        self.i_step = 0
+
         self.lattice_time = ""
         self.n_sites = 0
         self.bulk_counts = None
@@ -115,9 +136,19 @@ class SimController(object):
         self.vacuum_counts = self.sim.dope(self.empty_space_dopants, 'X')
 
     def run(self):
-        self.sim.run(self.kT, self.n_steps)
+        if type(self.kT) is list:
+            step_list = range(self.i_step, self.i_step + self.n_steps_per_frame)
+            exchanges_cumulative = [self.n_exchanges[0], self.n_exchanges[0] + self.n_exchanges[1]]
+            kT = np.interp(step_list, exchanges_cumulative, self.kT)
+            self.last_kT = np.average(kT)
+        else:
+            kT = self.kT
+            self.last_kT = kT
+        self.sim.run(kT, self.n_steps_per_frame)
+        self.i_step += self.n_steps_per_frame
 
     def write_frame(self, frames_path):
+        self.sim.lattice.info['T'] = self.last_kT
         self.sim.write(frames_path)
 
     def make_report_file(self, out_folder):
@@ -137,14 +168,26 @@ class SimController(object):
                                                 count=int(count),
                                                 perc=float(count)/sum_bulk*100)
 
-        vacuum_counts_string = ""
+        vacuum_counts_str = ""
         sum_vacuum = self.vacuum_counts[1].sum()
         for symbol, count in np.transpose(self.vacuum_counts):
-            vacuum_counts_string += elements_line.format(symbol=symbol,
+            vacuum_counts_str += elements_line.format(symbol=symbol,
                                                   count=int(count),
                                                   perc=float(count)/sum_vacuum*100)
 
-        lattice_counts = bulk_counts_header + bulk_counts_string + vacuum_counts_header + vacuum_counts_string
+        lattice_counts = bulk_counts_header + bulk_counts_string + vacuum_counts_header + vacuum_counts_str
+
+        if type(self.kT) is list:
+            kT_header_str = "kT [high, low]: "
+        else:
+            kT_header_str = "kT: "
+        kT_str = kT_header_str + str(self.kT)
+
+        if type(self.n_exchanges) is list:
+            exchanges_header_str = "Metropolis steps [high, fall, low]: "
+        else:
+            exchanges_header_str = "Metropolis steps: "
+        exchanges_str = exchanges_header_str + str(self.n_exchanges)
 
         report_filename = "sim_report.txt"
         report_path = os.path.join(out_folder, report_filename)
@@ -156,13 +199,11 @@ class SimController(object):
                                                   cell_x=cell_dimensions[0],
                                                   cell_y=cell_dimensions[1],
                                                   cell_z=cell_dimensions[2],
-
                                                   lattice_counts=lattice_counts,
-
                                                   particles_time=self.lattice_time,
 
-                                                  n_exchanges=self.n_exchanges,
-                                                  kT=self.kT,
+                                                  kT_str=kT_str,
+                                                  exchanges_str=exchanges_str,
                                                   n_frames=self.n_frames,
                                                   rand_ox=self.rand_ox,
                                                   sim_time=self.sim_time)
@@ -226,9 +267,12 @@ def make_config_file(conf_filename='default.in'):
     conf_writer['Job parameters'] = OrderedDict((('job_name', 'default'),
                                                  ('report_file', 'True')))
     conf_writer['Simulation parameters'] = OrderedDict((('sim_type', 'ase_sim'),
-                                                        ('kT', '300'),
-                                                        ('n_exchanges', '10000'),
-                                                        ('n_frames', '100'),
+                                                        ('kT_high', '500'),
+                                                        ('kT_low', '100'),
+                                                        ('n_steps_high', '10000'),
+                                                        ('n_steps_fall', '10000'),
+                                                        ('n_steps_low', '10000'),
+                                                        ('n_frames', '200'),
                                                         ('rand_ox', 'True')))
     conf_writer['Lattice parameters'] = OrderedDict((('bounds_scale', '10'),
                                                      ('x_scale', '1'),
@@ -353,20 +397,20 @@ for sim_ctrl in sim_controllers:
     start_time = time.time()
     print_time = time.time()
     # write first frame
-    idx = '1'.zfill(n_frames_width)
-    sim_frames_path = os.path.join(sim_frames_folder, sim_frames_filename.format(i=idx))
+    i_str = '1'.zfill(n_frames_width)
+    sim_frames_path = os.path.join(sim_frames_folder, sim_frames_filename.format(i=i_str))
     sim_ctrl.write_frame(sim_frames_path)
     for i in range(2, sim_ctrl.n_frames + 1):
         # simulate
         sim_ctrl.run()
         # write frame
-        idx = str(i).zfill(n_frames_width)
-        sim_frames_path = os.path.join(sim_frames_folder, sim_frames_filename.format(i=idx))
+        i_str = str(i).zfill(n_frames_width)
+        sim_frames_path = os.path.join(sim_frames_folder, sim_frames_filename.format(i=i_str))
         sim_ctrl.write_frame(sim_frames_path)
 
         # update console two times per second
         if time.time() - print_time >= 0.5:
-            print("\rWrote frame {i}/{n}".format(i=idx, n=str(sim_ctrl.n_frames)), end="")
+            print("\rWrote frame {i}/{n}".format(i=i_str, n=str(sim_ctrl.n_frames)), end="")
             print_time = time.time()
 
     print("\rWrote frame {}/{}".format(str(sim_ctrl.n_frames), str(sim_ctrl.n_frames)), flush=True)