init (NanoComp#955)

doc/docs/Python_User_Interface.md

@@ -912,7 +912,7 @@ where $G(t)$ is the current (not the dipole moment). In this formula, $\Delta f$
 
 ### CustomSource
 
-A user-specified source function $f(t)$. You can also specify start/end times at which point your current is set to zero whether or not your function is actually zero. These are optional, but you must specify an `end_time` explicitly if you want `run` functions like `until_after_sources` to work, since they need to know when your source turns off. For a demonstration of a [linear-chirped pulse](FAQ.md#how-do-i-create-a-chirped-pulse), see [`examples/chirped_pulse.py`](https://github.com/NanoComp/meep/blob/master/python/examples/chirped_pulse.py).
+A user-specified source function $f(t)$. You can also specify start/end times at which point your current is set to zero whether or not your function is actually zero. These are optional, but you must specify an `end_time` explicitly if you want `run` functions like `until_after_sources` to work, since they need to know when your source turns off. To use a custom source within an `EigenModeSource`, you must specify the `center_frequency` parameter, since Meep does not know the frequency content of the `CustomSource`. The resultant eigenmode is calculated at this frequency only. For a demonstration of a [linear-chirped pulse](FAQ.md#how-do-i-create-a-chirped-pulse), see [`examples/chirped_pulse.py`](https://github.com/NanoComp/meep/blob/master/python/examples/chirped_pulse.py).
 
 **`src_func` [`function`]**
 —
@@ -930,6 +930,10 @@ The end time for the source. Default is 10<sup>20</sup> (never turn off).
 —
 If `True`, the source is the integral of the current (the [dipole moment](https://en.wikipedia.org/wiki/Electric_dipole_moment)) which is guaranteed to be zero after the current turns off. In practice, there is little difference between integrated and non-integrated sources. Default is `False`.
 
+**`center_frequency` [`number`]**
+—
+Optional center frequency so that the `CustomSource` can be used within an `EigenModeSource`. Defaults to 0.
+
 <a name="fluxregion"></a>
 
 ### FluxRegion

python/source.py

@@ -71,12 +71,13 @@ def fourier_transform(self, freq):
 
 class CustomSource(SourceTime):
 
-    def __init__(self, src_func, start_time=-1.0e20, end_time=1.0e20, **kwargs):
+    def __init__(self, src_func, start_time=-1.0e20, end_time=1.0e20, center_frequency=0, **kwargs):
         super(CustomSource, self).__init__(**kwargs)
         self.src_func = src_func
         self.start_time = start_time
         self.end_time = end_time
-        self.swigobj = mp.custom_src_time(src_func, start_time, end_time)
+        self.center_frequency = center_frequency
+        self.swigobj = mp.custom_src_time(src_func, start_time, end_time, center_frequency)
         self.swigobj.is_integrated = self.is_integrated
 
 

python/tests/source.py

@@ -176,6 +176,53 @@ def test_amp_file_func(self):
         self.assertAlmostEqual(field_point_amp_file, field_point_amp_func, places=4)
         self.assertAlmostEqual(field_point_amp_arr, field_point_amp_func, places=4)
 
+class TestCustomEigenModeSource(unittest.TestCase):
+
+    def test_custom_em_source(self):
+        resolution = 20
+
+        dpml = 2
+        pml_layers = [mp.PML(thickness=dpml)]
+
+        sx = 40
+        sy = 12
+        cell_size = mp.Vector3(sx+2*dpml,sy)
+
+        v0 = 0.15  # pulse center frequency
+        a = 0.2*v0   # Gaussian envelope half-width
+        b = -0.1  # linear chirp rate (positive: up-chirp, negative: down-chirp)
+        t0 = 15   # peak time
+
+        chirp = lambda t: np.exp(1j*2*np.pi*v0*(t-t0)) * np.exp(-a*(t-t0)**2+1j*b*(t-t0)**2)
+
+        geometry = [mp.Block(center=mp.Vector3(0,0,0),size=mp.Vector3(mp.inf,1,mp.inf),material=mp.Medium(epsilon=12))]
+
+        kx = 0.4    # initial guess for wavevector in x-direction of eigenmode
+        kpoint = mp.Vector3(kx)
+        bnum = 1
+
+        sources = [mp.EigenModeSource(src=mp.CustomSource(src_func=chirp,center_frequency=v0),
+                            center=mp.Vector3(-0.5*sx + dpml + 1),
+                            size=mp.Vector3(y=sy),
+                            eig_kpoint=kpoint,
+                            eig_band=bnum,
+                            eig_parity=mp.EVEN_Y+mp.ODD_Z,
+                            eig_match_freq=True
+                            )]
+
+        sim = mp.Simulation(cell_size=cell_size,
+                            boundary_layers=pml_layers,
+                            resolution=resolution,
+                            k_point=mp.Vector3(),
+                            sources=sources,
+                            geometry=geometry,
+                            symmetries=[mp.Mirror(mp.Y)])
+
+        t = np.linspace(0,50,1e3)
+        sim.run(until=t0+50)
+
+        # For now, just check to make sure the simulation can run and the fields don't blow up.
+
 
 if __name__ == '__main__':
     unittest.main()

src/meep.hpp

@@ -918,8 +918,8 @@ class continuous_src_time : public src_time {
 class custom_src_time : public src_time {
 public:
   custom_src_time(std::complex<double> (*func)(double t, void *), void *data, double st = -infinity,
-                  double et = infinity)
-      : func(func), data(data), start_time(float(st)), end_time(float(et)) {}
+                  double et = infinity, std::complex<double> f = 0)
+      : func(func), data(data), start_time(float(st)), end_time(float(et)), freq(f) {}
   virtual ~custom_src_time() {}
 
   virtual std::complex<double> current(double time, double dt) const {
@@ -938,10 +938,13 @@ class custom_src_time : public src_time {
   virtual double last_time() const { return end_time; };
   virtual src_time *clone() const { return new custom_src_time(*this); }
   virtual bool is_equal(const src_time &t) const;
+  virtual std::complex<double> frequency() const { return freq; }
+  virtual void set_frequency(std::complex<double> f) { freq = f; }
 
 private:
   std::complex<double> (*func)(double t, void *);
   void *data;
+  std::complex<double> freq;
   double start_time, end_time;
 };