Add output functions and synchronized_magnetic (NanoComp#126)

python/examples/cyl-ellipsoid.py

@@ -26,16 +26,16 @@ def main():
 
     def print_stuff(sim_obj):
         v = mp.Vector3(4.13, 3.75, 0)
-        p = sim._get_field_point(src_cmpt, v)
-        print("t, Ez: {} {}+{}i".format(sim._round_time(), p.real, p.imag))
+        p = sim.get_field_point(src_cmpt, v)
+        print("t, Ez: {} {}+{}i".format(sim.round_time(), p.real, p.imag))
 
     sim.run(mp.at_beginning(mp.output_epsilon),
             mp.at_every(0.25, print_stuff),
             mp.at_end(print_stuff),
             mp.at_end(mp.output_efield_z),
             until=23)
 
-    print("stopped at meep time = {}".format(sim._round_time()))
+    print("stopped at meep time = {}".format(sim.round_time()))
 
 
 if __name__ == '__main__':

python/meep.i

@@ -403,10 +403,21 @@ extern boolean point_in_objectp(vector3 p, GEOMETRIC_OBJECT o);
         in_volume,
         interpolate,
         output_epsilon,
+        output_mu,
+        output_hpwr,
+        output_dpwr,
+        output_tot_pwr,
         output_hfield_z,
         output_efield_z,
+        output_poynting,
+        output_poynting_x,
+        output_poynting_y,
+        output_poynting_z,
+        output_poynting_r,
+        output_poynting_p,
         py_v3_to_vec,
         stop_when_fields_decayed,
+        synchronized_magnetic,
         to_appended
     )
     from .source import (

python/simulation.py

@@ -185,7 +185,7 @@ def _collect1(c, pt):
             def _collect2(sim):
                 self.data_dt = sim.meep_time() - self.t0
                 self.t0 = sim.meep_time()
-                self.data.append(sim._get_field_point(c, pt))
+                self.data.append(sim.get_field_point(c, pt))
             return _collect2
         return _collect1
 
@@ -376,13 +376,13 @@ def meep_time(self):
             self._init_fields()
         return self.fields.time()
 
-    def _round_time(self):
+    def round_time(self):
         if self.fields is None:
             self._init_fields()
 
         return self.fields.round_time()
 
-    def _get_field_point(self, c, pt):
+    def get_field_point(self, c, pt):
         v3 = py_v3_to_vec(self.dimensions, pt, self.is_cylindrical)
         return self.fields.get_field_from_comp(c, v3)
 
@@ -401,10 +401,10 @@ def _run_until(self, cond, step_funcs):
 
         if isinstance(cond, numbers.Number):
             stop_time = cond
-            t0 = self._round_time()
+            t0 = self.round_time()
 
             def stop_cond(sim):
-                return sim._round_time() >= t0 + stop_time
+                return sim.round_time() >= t0 + stop_time
 
             cond = stop_cond
 
@@ -435,10 +435,10 @@ def _run_sources_until(self, cond, step_funcs):
         ts = self.fields.last_source_time()
 
         if isinstance(cond, numbers.Number):
-            new_cond = (ts - self._round_time()) + cond
+            new_cond = (ts - self.round_time()) + cond
         else:
             def f(sim):
-                return cond(sim) and sim._round_time() >= ts
+                return cond(sim) and sim.round_time() >= ts
             new_cond = f
 
         self._run_until(new_cond, step_funcs)
@@ -632,6 +632,26 @@ def output_component(self, c, h5file=None):
                 self.last_eps_filename = nm
             self.output_h5_hook(nm)
 
+    def output_components(self, fname, *components):
+        if self.fields is None:
+            raise RuntimeError("Fields must be initialized before calling output_component")
+
+        if self.output_append_h5 is None:
+            f = self.fields.open_h5file(fname, mp.h5file.WRITE, self._get_filename_prefix(), True)
+        else:
+            f = None
+
+        for c in components:
+            self.output_component(c, h5file=f)
+            if self.output_append_h5 is None:
+                f.prevent_deadlock()
+
+        if self.output_append_h5 is None:
+            f = None
+
+        if self.output_append_h5 is None:
+            self.output_h5_hook(self.fields.h5file_name(fname, self._get_filename_prefix(), True))
+
     def change_k_point(self, k):
         self.k_point = k
 
@@ -753,15 +773,15 @@ def _true(sim, todo):
 def after_sources(*step_funcs):
     def _after_sources(sim, todo):
         time = sim.fields.last_source_time()
-        if sim._round_time() >= time:
+        if sim.round_time() >= time:
             for func in step_funcs:
                 _eval_step_func(sim, func, todo)
     return _after_sources
 
 
 def after_time(t, *step_funcs):
     def _after_t(sim):
-        return sim._round_time() >= t
+        return sim.round_time() >= t
     return _when_true_funcs(_after_t, *step_funcs)
 
 
@@ -790,7 +810,7 @@ def at_every(dt, *step_funcs):
     closure = {'tlast': 0.0}
 
     def _every(sim, todo):
-        t = sim._round_time()
+        t = sim.round_time()
         if todo == 'finish' or t >= closure['tlast'] + dt + (-0.5 * sim.fields.dt):
             for func in step_funcs:
                 _eval_step_func(sim, func, todo)
@@ -800,7 +820,7 @@ def _every(sim, todo):
 
 def before_time(t, *step_funcs):
     def _before_t(sim):
-        return sim._round_time() < t
+        return sim.round_time() < t
     return _when_true_funcs(_before_t, *step_funcs)
 
 
@@ -809,7 +829,7 @@ def during_sources(*step_funcs):
 
     def _during_sources(sim, todo):
         time = sim.fields.last_source_time()
-        if sim._round_time() < time:
+        if sim.round_time() < time:
             for func in step_funcs:
                 _eval_step_func(sim, func, 'step')
         elif closure['finished'] is False:
@@ -871,15 +891,15 @@ def stop_when_fields_decayed(dt, c, pt, decay_by):
     }
 
     def _stop(sim):
-        fabs = abs(sim._get_field_point(c, pt)) * abs(sim._get_field_point(c, pt))
+        fabs = abs(sim.get_field_point(c, pt)) * abs(sim.get_field_point(c, pt))
         closure['cur_max'] = max(closure['cur_max'], fabs)
 
-        if sim._round_time() <= dt + closure['t0']:
+        if sim.round_time() <= dt + closure['t0']:
             return False
         else:
             old_cur = closure['cur_max']
             closure['cur_max'] = 0
-            closure['t0'] = sim._round_time()
+            closure['t0'] = sim.round_time()
             closure['max_abs'] = max(closure['max_abs'], old_cur)
             if closure['max_abs'] != 0:
                 fmt = "field decay(t = {}): {} / {} = {}"
@@ -888,6 +908,15 @@ def _stop(sim):
     return _stop
 
 
+def synchronized_magnetic(*step_funcs):
+    def _sync(sim, todo):
+        sim.fields.synchronize_magnetic_fields()
+        for f in step_funcs:
+            _eval_step_func(sim, f, todo)
+        sim.fields.restore_magnetic_fields()
+    return _sync
+
+
 def display_csv(sim, name, data):
     for d in data:
         display_run_data(sim, name, d)
@@ -931,6 +960,22 @@ def output_epsilon(sim):
     sim.output_component(mp.Dielectric)
 
 
+def output_mu(sim):
+    sim.output_component(mp.Permeability)
+
+
+def output_hpwr(sim):
+    sim.output_component(mp.H_EnergyDensity)
+
+
+def output_dpwr(sim):
+    sim.output_component(mp.D_EnergyDensity)
+
+
+def output_tot_pwr(sim):
+    sim.output_component(mp.EnergyDensity)
+
+
 # TODO(chogan): Write rest of convenience functions generated in meep.scm:define-output-field
 def output_hfield_z(sim):
     sim.output_component(mp.Hz)
@@ -940,6 +985,30 @@ def output_efield_z(sim):
     sim.output_component(mp.Ez)
 
 
+def output_poynting(sim):
+    sim.output_components('s', mp.Sx, mp.Sy, mp.Sz, mp.Sr, mp.Sp)
+
+
+def output_poynting_x(sim):
+    sim.output_component(mp.Sx)
+
+
+def output_poynting_y(sim):
+    sim.output_component(mp.Sy)
+
+
+def output_poynting_z(sim):
+    sim.output_component(mp.Sz)
+
+
+def output_poynting_r(sim):
+    sim.output_component(mp.Sr)
+
+
+def output_poynting_p(sim):
+    sim.output_component(mp.Sp)
+
+
 def get_flux_freqs(f):
     return np.linspace(f.freq_min, f.freq_min + f.dfreq * f.Nfreq, num=f.Nfreq, endpoint=False).tolist()
 

python/tests/cyl_ellipsoid.py

@@ -37,8 +37,8 @@ def init(self):
 
         def print_stuff(sim_obj):
             v = mp.Vector3(4.13, 3.75, 0)
-            p = self.sim._get_field_point(self.src_cmpt, v)
-            print("t, Ez: {} {}+{}i".format(self.sim._round_time(), p.real, p.imag))
+            p = self.sim.get_field_point(self.src_cmpt, v)
+            print("t, Ez: {} {}+{}i".format(self.sim.round_time(), p.real, p.imag))
         self.print_stuff = print_stuff
 
     def run_simulation(self):

python/tests/ring.py

@@ -56,7 +56,7 @@ def test_harminv(self):
         self.assertAlmostEqual(m1.amp.real, -0.00304951667301, places=4)
         self.assertAlmostEqual(m1.amp.imag, -0.00153192946717, places=4)
 
-        fp = self.sim._get_field_point(mp.Ez, mp.Vector3(1, 1))
+        fp = self.sim.get_field_point(mp.Ez, mp.Vector3(1, 1))
         self.assertAlmostEqual(fp, -0.08185972142450348)
 
 if __name__ == '__main__':