Bugfix model.wavenumber for ab=36 (0s)

CHANGELOG.rst

@@ -16,32 +16,33 @@ Changelog
 - Unified DLF arguments
   [`empymod#10 <https://github.com/empymod/empymod/issues/10>`_].
 
-    These changes are backwards compatible for all main modelling routines in
-    ``empymod.model``. However, they are not backwards compatible for the
-    following routines:
-    - ``empymod.model.fem`` (removed ``use_spline``),
-    - ``empymod.transform.fht`` (removed ``use_spline``),
-    - ``empymod.transform.hqwe`` (removed ``use_spline``),
-    - ``empymod.transform.quad`` (removed ``use_spline``),
-    - ``empymod.transform.dlf`` (``lagged``, ``splined`` repl. by
-      ``pts_per_dec``),
-    - ``empymod.utils.check_opt`` (no longer returns ``use_spline``),
-    - ``empymod.utils.check_hankel`` (changes in ``pts_per_dec``), and
-    - ``empymod.utils.check_time`` (changes in ``pts_per_dec``).
-
-    The function ``empymod.utils.spline_backwards_hankel`` can be used for
-    backwards compatibility.
-
-    Now the Hankel and Fourier DLF have the same behaviour for
-    ``pts_per_dec``:
-    - ``pts_per_dec = 0``: Standard DLF,
-    - ``pts_per_dec < 0``: Lagged Convolution DLF, and
-    - ``pts_per_dec > 0``: Splined DLF.
-
-    **There is one exception** which is not backwards compatible: Before,
-    if ``opt=None`` and ``htarg={pts_per_dec: != 0}``, the ``pts_per_dec`` was
-    not used for the FHT and the QWE. New, this will be used according to the
-    above definitions.
+  These changes are backwards compatible for all main modelling routines in
+  ``empymod.model``. However, they are not backwards compatible for the
+  following routines:
+  - ``empymod.model.fem`` (removed ``use_spline``),
+  - ``empymod.transform.fht`` (removed ``use_spline``),
+  - ``empymod.transform.hqwe`` (removed ``use_spline``),
+  - ``empymod.transform.quad`` (removed ``use_spline``),
+  - ``empymod.transform.dlf`` (``lagged``, ``splined`` repl. by
+    ``pts_per_dec``),
+  - ``empymod.utils.check_opt`` (no longer returns ``use_spline``),
+  - ``empymod.utils.check_hankel`` (changes in ``pts_per_dec``), and
+  - ``empymod.utils.check_time`` (changes in ``pts_per_dec``).
+
+  The function ``empymod.utils.spline_backwards_hankel`` can be used for
+  backwards compatibility.
+
+  Now the Hankel and Fourier DLF have the same behaviour for ``pts_per_dec``:
+  - ``pts_per_dec = 0``: Standard DLF,
+  - ``pts_per_dec < 0``: Lagged Convolution DLF, and
+  - ``pts_per_dec > 0``: Splined DLF.
+
+  **There is one exception** which is not backwards compatible: Before, if
+  ``opt=None`` and ``htarg={pts_per_dec: != 0}``, the ``pts_per_dec`` was not
+  used for the FHT and the QWE. New, this will be used according to the above
+  definitions.
+
+- Bugfix in ``model.wavenumber`` for ``ab=36`` (zeroes).
 
 
 v1.5.2 - *2018-04-25*

empymod/model.py

@@ -1389,7 +1389,7 @@ def wavenumber(src, rec, depth, res, freq, wavenumber, ab=11, aniso=None,
 
     # Check frequency => get etaH, etaV, zetaH, and zetaV
     f = check_frequency(freq, res, aniso, epermH, epermV, mpermH, mpermV, verb)
-    _, etaH, etaV, zetaH, zetaV = f  # (output freq not required)
+    freq, etaH, etaV, zetaH, zetaV = f
 
     # Check src-rec configuration
     # => Get flags if src or rec or both are magnetic (msrc, mrec)
@@ -1409,6 +1409,11 @@ def wavenumber(src, rec, depth, res, freq, wavenumber, ab=11, aniso=None,
 
     # === 3. EM-FIELD CALCULATION ============
 
+    # If <ab> = 36 (or 63), field is zero
+    if ab_calc in [36, ]:
+        out = np.zeros((freq.size, off.size, wavenumber.size), dtype=complex)
+        return np.squeeze(out), np.squeeze(out)
+
     # Calculate wavenumber response
     PJ0, PJ1, PJ0b = kernel.wavenumber(zsrc, zrec, lsrc, lrec, depth, etaH,
                                        etaV, zetaH, zetaV,

tests/test_model.py

@@ -555,6 +555,12 @@ def test_wavenumber():
     assert_allclose(w_res0, res['PJ0'])
     assert_allclose(w_res1, res['PJ1'])
 
+    # Check that ab=36 returns zeros
+    res['inp']['ab'] = 36
+    w_res0, w_res1 = wavenumber(**res['inp'])
+    assert_allclose(w_res0, np.zeros(res['PJ0'].shape, dtype=complex))
+    assert_allclose(w_res1, np.zeros(res['PJ1'].shape, dtype=complex))
+
 
 def test_fem():
     # Just ensure functionality stays the same, with one example.