New GaussianPulse.from_frequency_range for maximizing amplitude in frequency range

Author: weiliangjin2021

CHANGELOG.md

@@ -12,6 +12,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Periodic repetition of EME subgrids via `num_reps` or `EMEPeriodicitySweep`.
 - Methods `EMEExplicitGrid.from_structures` and `EMECompositeGrid.from_structure_groups` to place EME cell boundaries at structure bounds.
 - 'ModeSimulation' now supports 'PermittivityMonitor'.
+- Classmethod `from_frequency_range` in `GaussianPulse` for generating a pulse whose amplitude in the frequency_range [fmin, fmax] is maximized, which is particularly useful for running broadband simulations.
 
 ### Changed
 - Performance enhancement for adjoint gradient calculations by optimizing field interpolation.

tests/test_components/test_source.py

@@ -93,6 +93,47 @@ def test_source_times():
     assert abs(dc_comp) ** 2 > 1e-32
 
 
+def test_gaussian_from_frequency_range():
+    # error with negative fmin
+    with pytest.raises(ValueError):
+        _ = td.GaussianPulse.from_frequency_range(fmin=-1e10, fmax=1e10)
+    # error with fmin = 0
+    with pytest.raises(ValueError):
+        _ = td.GaussianPulse.from_frequency_range(fmin=0, fmax=1e10)
+    # error with fmin >= fmax
+    with pytest.raises(ValueError):
+        _ = td.GaussianPulse.from_frequency_range(fmin=1e10, fmax=0.9e10)
+
+    fmin = 1e9
+    fmax = 20e9
+    # dc component on
+    g = td.GaussianPulse.from_frequency_range(fmin=fmin, fmax=fmax, remove_dc_component=False)
+    assert g.freq0 == 0.5 * (fmin + fmax)
+    assert g.fwidth == 0.5 * (fmax - fmin)
+
+    # dc component removed
+    g1 = td.GaussianPulse.from_frequency_range(fmin=fmin, fmax=fmax, remove_dc_component=True)
+    # default to dc removed
+    g2 = td.GaussianPulse.from_frequency_range(fmin=fmin, fmax=fmax)
+    assert g2.remove_dc_component
+
+    # 1) broadband: assert enough amplitude at fmin and fmax
+    time = np.linspace(0, 5 / fmin, 10001)
+    freqs = np.linspace(fmin, fmax, 101)
+    spectrum = np.abs(g2.spectrum(time, freqs, time[1] - time[0]))
+    max_amp = np.max(spectrum)
+    assert spectrum[0] / max_amp > 0.1
+    assert spectrum[-1] / max_amp > 0.1
+
+    # 2) narrow band: close to regular Gaussian result
+    fmin = 10e9
+    bandwidth = 1e6
+    fmax = fmin + 2 * bandwidth
+    g = td.GaussianPulse.from_frequency_range(fmin=fmin, fmax=fmax)
+    assert abs(g.fwidth - bandwidth) / bandwidth < 1e-4
+    assert abs(g.freq0 - fmin) / fmin < 1e-4
+
+
 def test_dipole():
     g = td.GaussianPulse(freq0=1e12, fwidth=0.1e12)
     _ = td.PointDipole(center=(1, 2, 3), source_time=g, polarization="Ex", interpolate=True)

tidy3d/components/source/time.py

@@ -191,6 +191,48 @@ def from_amp_complex(cls, amp: complex, **kwargs) -> GaussianPulse:
         phase = np.angle(amp)
         return cls(amplitude=amplitude, phase=phase, **kwargs)
 
+    @classmethod
+    def from_frequency_range(
+        cls, fmin: pydantic.PositiveFloat, fmax: pydantic.PositiveFloat, **kwargs
+    ) -> GaussianPulse:
+        """Create a ``GaussianPulse`` that maximizes its amplitude in the frequency range [fmin, fmax].
+
+        Parameters
+        ----------
+        fmin : float
+            Lower bound of frequency of interest.
+        fmax : float
+            Upper bound of frequency of interest.
+        kwargs : dict
+            Keyword arguments passed to ``GaussianPulse()``, excluding ``freq0`` & ``fwidth``.
+
+        Returns
+        -------
+        GaussianPulse
+            A ``GaussianPulse`` that maximizes its amplitude in the frequency range [fmin, fmax].
+        """
+        # validate that fmin and fmax must positive, and fmax > fmin
+        if fmin <= 0:
+            raise ValidationError("'fmin' must be positive.")
+        if fmax <= fmin:
+            raise ValidationError("'fmax' must be greater than 'fmin'.")
+
+        # frequency range and center
+        freq_range = fmax - fmin
+        freq_center = (fmax + fmin) / 2.0
+
+        # If remove_dc_component=False, simply return the standard GaussianPulse parameters
+        if kwargs.get("remove_dc_component", True) is False:
+            return cls(freq0=freq_center, fwidth=freq_range / 2.0, **kwargs)
+
+        # If remove_dc_component=True, the Gaussian pulse is distorted
+        kwargs.update({"remove_dc_component": True})
+        log_ratio = np.log(fmax / fmin)
+        coeff = ((1 + log_ratio**2) ** 0.5 - 1) / 2.0
+        freq0 = freq_center - coeff / log_ratio * freq_range
+        fwidth = freq_range / log_ratio * coeff**0.5
+        return cls(freq0=freq0, fwidth=fwidth, **kwargs)
+
 
 class ContinuousWave(Pulse):
     """Source time dependence that ramps up to continuous oscillation

tidy3d/plugins/smatrix/component_modelers/terminal.py

@@ -31,7 +31,7 @@
 from ..ports.coaxial_lumped import CoaxialLumpedPort
 from ..ports.rectangular_lumped import LumpedPort
 from ..ports.wave import WavePort
-from .base import FWIDTH_FRAC, AbstractComponentModeler, TerminalPortType
+from .base import AbstractComponentModeler, TerminalPortType
 
 
 class TerminalComponentModeler(AbstractComponentModeler):
@@ -164,11 +164,8 @@ def sim_dict(self) -> Dict[str, Simulation]:
     @cached_property
     def _source_time(self):
         """Helper to create a time domain pulse for the frequency range of interest."""
-        freq0 = np.mean(self.freqs)
-        fdiff = max(self.freqs) - min(self.freqs)
-        fwidth = max(fdiff, freq0 * FWIDTH_FRAC)
-        return GaussianPulse(
-            freq0=freq0, fwidth=fwidth, remove_dc_component=self.remove_dc_component
+        return GaussianPulse.from_frequency_range(
+            fmin=min(self.freqs), fmax=max(self.freqs), remove_dc_component=self.remove_dc_component
         )
 
     def _construct_smatrix(self) -> TerminalPortDataArray: