diff --git a/.github/images/anim_fdtd_metadiffuser.gif b/.github/images/anim_fdtd_metadiffuser.gif
new file mode 100644
index 000000000..846cf0edd
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser.gif differ
diff --git a/.github/images/anim_fdtd_metadiffuser.webm b/.github/images/anim_fdtd_metadiffuser.webm
new file mode 100644
index 000000000..020a0b67b
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser.webm differ
diff --git a/.github/images/anim_fdtd_metadiffuser_dark.gif b/.github/images/anim_fdtd_metadiffuser_dark.gif
new file mode 100644
index 000000000..123b8b65a
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser_dark.gif differ
diff --git a/.github/images/anim_fdtd_metadiffuser_dark.webm b/.github/images/anim_fdtd_metadiffuser_dark.webm
new file mode 100644
index 000000000..a2fb4559c
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser_dark.webm differ
diff --git a/.github/images/anim_fdtd_metadiffuser_dark_poster.jpg b/.github/images/anim_fdtd_metadiffuser_dark_poster.jpg
new file mode 100644
index 000000000..5561e6194
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser_dark_poster.jpg differ
diff --git a/.github/images/anim_fdtd_metadiffuser_es.webm b/.github/images/anim_fdtd_metadiffuser_es.webm
new file mode 100644
index 000000000..f5b5bcf59
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser_es.webm differ
diff --git a/.github/images/anim_fdtd_metadiffuser_es_dark.webm b/.github/images/anim_fdtd_metadiffuser_es_dark.webm
new file mode 100644
index 000000000..6ea66dc6e
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser_es_dark.webm differ
diff --git a/.github/images/anim_fdtd_metadiffuser_es_dark_poster.jpg b/.github/images/anim_fdtd_metadiffuser_es_dark_poster.jpg
new file mode 100644
index 000000000..11f3e8518
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser_es_dark_poster.jpg differ
diff --git a/.github/images/anim_fdtd_metadiffuser_es_poster.jpg b/.github/images/anim_fdtd_metadiffuser_es_poster.jpg
new file mode 100644
index 000000000..2197537b6
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser_es_poster.jpg differ
diff --git a/.github/images/anim_fdtd_metadiffuser_poster.jpg b/.github/images/anim_fdtd_metadiffuser_poster.jpg
new file mode 100644
index 000000000..d74cd2bf4
Binary files /dev/null and b/.github/images/anim_fdtd_metadiffuser_poster.jpg differ
diff --git a/.github/images/metadiffuser_geometry.svg b/.github/images/metadiffuser_geometry.svg
new file mode 100644
index 000000000..a61cb76a1
--- /dev/null
+++ b/.github/images/metadiffuser_geometry.svg
@@ -0,0 +1,425 @@
+
+
+
+
+
+
+ image/svg+xml
+
+
+ Matplotlib v3.11.1, https://matplotlib.org/
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ 1
+
+
+ 2
+
+
+ 3
+
+
+ 4
+
+
+ 5
+
+
+ Incident sound
+
+
+ Metadiffuser cross-section (one period)
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ 350 mm
+
+
+ 20 mm
+
+
+ 70 mm
+
+
+ 14.7 mm
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/.github/images/metadiffuser_geometry_dark.svg b/.github/images/metadiffuser_geometry_dark.svg
new file mode 100644
index 000000000..d5ed2c141
--- /dev/null
+++ b/.github/images/metadiffuser_geometry_dark.svg
@@ -0,0 +1,425 @@
+
+
+
+
+
+
+ image/svg+xml
+
+
+ Matplotlib v3.11.1, https://matplotlib.org/
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ 1
+
+
+ 2
+
+
+ 3
+
+
+ 4
+
+
+ 5
+
+
+ Incident sound
+
+
+ Metadiffuser cross-section (one period)
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ 350 mm
+
+
+ 20 mm
+
+
+ 70 mm
+
+
+ 14.7 mm
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/.github/images/metadiffuser_geometry_es.svg b/.github/images/metadiffuser_geometry_es.svg
new file mode 100644
index 000000000..2b1846ae1
--- /dev/null
+++ b/.github/images/metadiffuser_geometry_es.svg
@@ -0,0 +1,425 @@
+
+
+
+
+
+
+ image/svg+xml
+
+
+ Matplotlib v3.11.1, https://matplotlib.org/
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ 1
+
+
+ 2
+
+
+ 3
+
+
+ 4
+
+
+ 5
+
+
+ Sonido incidente
+
+
+ Sección del metadifusor (un periodo)
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ 350 mm
+
+
+ 20 mm
+
+
+ 70 mm
+
+
+ 14,7 mm
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/.github/images/metadiffuser_geometry_es_dark.svg b/.github/images/metadiffuser_geometry_es_dark.svg
new file mode 100644
index 000000000..358704cc2
--- /dev/null
+++ b/.github/images/metadiffuser_geometry_es_dark.svg
@@ -0,0 +1,425 @@
+
+
+
+
+
+
+ image/svg+xml
+
+
+ Matplotlib v3.11.1, https://matplotlib.org/
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ 1
+
+
+ 2
+
+
+ 3
+
+
+ 4
+
+
+ 5
+
+
+ Sonido incidente
+
+
+ Sección del metadifusor (un periodo)
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ 350 mm
+
+
+ 20 mm
+
+
+ 70 mm
+
+
+ 14,7 mm
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/.github/images/metadiffuser_polar.svg b/.github/images/metadiffuser_polar.svg
new file mode 100644
index 000000000..737e0feea
--- /dev/null
+++ b/.github/images/metadiffuser_polar.svg
@@ -0,0 +1,457 @@
+
+
+
+
+
+
+ image/svg+xml
+
+
+ Matplotlib v3.11.1, https://matplotlib.org/
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ -90°
+
+
+
+
+
+
+ -60°
+
+
+
+
+
+
+ -30°
+
+
+
+
+
+
+ 0°
+
+
+
+
+
+
+ 30°
+
+
+
+
+
+
+ 60°
+
+
+
+
+
+
+ 90°
+
+
+
+
+
+
+
+
+ −35
+
+
+
+
+
+
+ −30
+
+
+
+
+
+
+ −25
+
+
+
+
+
+
+ −20
+
+
+
+
+
+
+ −15
+
+
+
+
+
+
+ −10
+
+
+
+
+
+
+ −5
+
+
+
+
+
+
+ 0
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ The 2 cm metadiffuser scatters like the 27 cm QRD (2 kHz)
+
+
+
+
+
+
+
+ Metadiffuser, panel 2 cm
+
+
+
+
+ QRD, wells up to 27.4 cm
+
+
+
+
+
+
+
+
+
diff --git a/.github/images/metadiffuser_polar_dark.svg b/.github/images/metadiffuser_polar_dark.svg
new file mode 100644
index 000000000..c85421bc1
--- /dev/null
+++ b/.github/images/metadiffuser_polar_dark.svg
@@ -0,0 +1,457 @@
+
+
+
+
+
+
+ image/svg+xml
+
+
+ Matplotlib v3.11.1, https://matplotlib.org/
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ -90°
+
+
+
+
+
+
+ -60°
+
+
+
+
+
+
+ -30°
+
+
+
+
+
+
+ 0°
+
+
+
+
+
+
+ 30°
+
+
+
+
+
+
+ 60°
+
+
+
+
+
+
+ 90°
+
+
+
+
+
+
+
+
+ −35
+
+
+
+
+
+
+ −30
+
+
+
+
+
+
+ −25
+
+
+
+
+
+
+ −20
+
+
+
+
+
+
+ −15
+
+
+
+
+
+
+ −10
+
+
+
+
+
+
+ −5
+
+
+
+
+
+
+ 0
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ The 2 cm metadiffuser scatters like the 27 cm QRD (2 kHz)
+
+
+
+
+
+
+
+ Metadiffuser, panel 2 cm
+
+
+
+
+ QRD, wells up to 27.4 cm
+
+
+
+
+
+
+
+
+
diff --git a/.github/images/metadiffuser_polar_es.svg b/.github/images/metadiffuser_polar_es.svg
new file mode 100644
index 000000000..73cb6a2f5
--- /dev/null
+++ b/.github/images/metadiffuser_polar_es.svg
@@ -0,0 +1,457 @@
+
+
+
+
+
+
+ image/svg+xml
+
+
+ Matplotlib v3.11.1, https://matplotlib.org/
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ -90°
+
+
+
+
+
+
+ -60°
+
+
+
+
+
+
+ -30°
+
+
+
+
+
+
+ 0°
+
+
+
+
+
+
+ 30°
+
+
+
+
+
+
+ 60°
+
+
+
+
+
+
+ 90°
+
+
+
+
+
+
+
+
+ −35
+
+
+
+
+
+
+ −30
+
+
+
+
+
+
+ −25
+
+
+
+
+
+
+ −20
+
+
+
+
+
+
+ −15
+
+
+
+
+
+
+ −10
+
+
+
+
+
+
+ −5
+
+
+
+
+
+
+ 0
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ El metadifusor de 2 cm dispersa como el QRD de 27 cm (2 kHz)
+
+
+
+
+
+
+
+ Metadifusor, panel de 2 cm
+
+
+
+
+ QRD, pozos de hasta 27,4 cm
+
+
+
+
+
+
+
+
+
diff --git a/.github/images/metadiffuser_polar_es_dark.svg b/.github/images/metadiffuser_polar_es_dark.svg
new file mode 100644
index 000000000..3f2ed8111
--- /dev/null
+++ b/.github/images/metadiffuser_polar_es_dark.svg
@@ -0,0 +1,457 @@
+
+
+
+
+
+
+ image/svg+xml
+
+
+ Matplotlib v3.11.1, https://matplotlib.org/
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ -90°
+
+
+
+
+
+
+ -60°
+
+
+
+
+
+
+ -30°
+
+
+
+
+
+
+ 0°
+
+
+
+
+
+
+ 30°
+
+
+
+
+
+
+ 60°
+
+
+
+
+
+
+ 90°
+
+
+
+
+
+
+
+
+ −35
+
+
+
+
+
+
+ −30
+
+
+
+
+
+
+ −25
+
+
+
+
+
+
+ −20
+
+
+
+
+
+
+ −15
+
+
+
+
+
+
+ −10
+
+
+
+
+
+
+ −5
+
+
+
+
+
+
+ 0
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+ El metadifusor de 2 cm dispersa como el QRD de 27 cm (2 kHz)
+
+
+
+
+
+
+
+ Metadifusor, panel de 2 cm
+
+
+
+
+ QRD, pozos de hasta 27,4 cm
+
+
+
+
+
+
+
+
+
diff --git a/docs/api-reference.md b/docs/api-reference.md
index bac4df137..085c769cb 100644
--- a/docs/api-reference.md
+++ b/docs/api-reference.md
@@ -833,6 +833,12 @@ cycle and warn on use; they are removed in 4.0.
| `quadratic_residue_sequence` | `function` | **Quadratic residue sequence s_n (Cox & D'Antonio Eq. 10.2).**
+Show the code for this figure
+
+```python
+import numpy as np
+from phonometry import (
+ HelmholtzResonator,
+ MetadiffuserWell,
+ materials,
+ metadiffuser_polar_response,
+)
+
+mm = 1e-3
+rows = [
+ (14.7, 13.0, 16.4, 6.2, 9.0),
+ (30.9, 9.1, 4.3, 3.5, 9.0),
+ (30.9, 9.1, 4.3, 3.5, 9.0),
+ (15.7, 13.3, 17.0, 6.3, 9.0),
+ (20.3, 18.0, 20.7, 3.2, 9.0),
+]
+wells = [
+ MetadiffuserWell(
+ h * mm,
+ 2 * (HelmholtzResonator(ln * mm, wn * mm, lc * mm, wc * mm),),
+ )
+ for h, ln, lc, wn, wc in rows
+]
+
+# The metadiffuser panel and the QRD it was tuned to at 2 kHz, both with
+# six repetitions of the period.
+meta = metadiffuser_polar_response(2000.0, wells, depth=0.02, period=0.07,
+ periods=6)
+sequence = np.roll(materials.quadratic_residue_sequence(5), -1)
+depths = sequence * (343.0 / 500.0) / (2 * 5)
+qrd = materials.predict_diffuser_polar_response(
+ 0.07, 2000.0, depths=depths, periods=6, include_obliquity=False,
+)
+
+ax = meta.plot(marker="", linewidth=2.2, label="Metadiffuser, panel 2 cm")
+qrd.plot(ax=ax, marker="", linewidth=1.6, linestyle="--",
+ label="QRD, wells up to 27.4 cm")
+ax.legend(loc="lower center")
+```
+
+
+
+The far-field overlay above is the frequency-domain summary; the FDTD
+animation below meshes both panels for real (the metadiffuser at 0.25 mm, slits,
+necks and cavities included) and lets the same 2 kHz wavefront hit them: the
+27 cm QRD and the 2 cm panel throw out nearly the same scattered fan.
+
+ `` and streams raw
+ frames on stdin; the shim forwards stdin over ssh into the NVENC
+ container, which writes the WebM to a file in the remote work
+ directory (the matroska muxer needs a seekable target), and the shim
+ then fetches it back with scp and cleans up.
+ """
+ import shlex
+ import tempfile
+
+ q_target = shlex.quote(target)
+ q_image = shlex.quote(image)
+ q_dir = shlex.quote(workdir)
+ script = f"""#!/bin/bash
+set -euo pipefail
+out="${{@: -1}}"
+args=("${{@:1:$#-1}}")
+rid="enc-$$-$RANDOM.webm"
+quoted=$(printf ' %q' "${{args[@]}}")
+ssh -o BatchMode=yes -- {q_target} \
+ "mkdir -p {q_dir} && docker run --rm -i --gpus all \
+ -v {q_dir}:/work {q_image}$quoted /work/$rid"
+scp -q -- {q_target}:{q_dir}/"$rid" "$out"
+ssh -o BatchMode=yes -- {q_target} "rm -f {q_dir}/$rid"
+"""
+ with tempfile.NamedTemporaryFile(
+ "w", suffix=".sh", prefix="ffmpeg-gpu-", delete=False) as handle:
+ handle.write(script)
+ os.chmod(handle.name, 0o755)
+ return handle.name
+
+
def _save_animation(anim: Any, fig: Any, output_dir: str, stem: str,
make_gif: bool = True, *, fps: int | None = None,
gif_fps: int | None = None) -> None:
@@ -12573,13 +12761,30 @@ def _save_animation(anim: Any, fig: Any, output_dir: str, stem: str,
from matplotlib.animation import FFMpegWriter
webm = _anim_path(output_dir, stem, "webm")
- writer = FFMpegWriter(
- fps=_ANIM_FPS if fps is None else fps, codec="libvpx-vp9",
- extra_args=_vp9_extra_args(),
- )
- with plt.rc_context({"savefig.bbox": "standard"}):
- anim.save(webm, writer=writer, dpi=_ANIM_DPI,
- savefig_kwargs={"facecolor": fig.get_facecolor()})
+ rate = _ANIM_FPS if fps is None else fps
+ gpu = _gpu_encode_target()
+ saved = False
+ if gpu is not None:
+ wrapper = _write_gpu_ffmpeg_wrapper(gpu["target"], gpu["image"],
+ gpu["workdir"])
+ try:
+ writer = FFMpegWriter(fps=rate, codec="av1_nvenc",
+ extra_args=_av1_nvenc_extra_args())
+ with plt.rc_context({"savefig.bbox": "standard",
+ "animation.ffmpeg_path": wrapper}):
+ anim.save(webm, writer=writer, dpi=_ANIM_DPI,
+ savefig_kwargs={"facecolor": fig.get_facecolor()})
+ saved = True
+ except (RuntimeError, subprocess.CalledProcessError, OSError) as exc:
+ print(f" [gpu-encode] {stem}: {exc}; falling back to VP9")
+ finally:
+ os.remove(wrapper)
+ if not saved:
+ writer = FFMpegWriter(fps=rate, codec="libvpx-vp9",
+ extra_args=_vp9_extra_args())
+ with plt.rc_context({"savefig.bbox": "standard"}):
+ anim.save(webm, writer=writer, dpi=_ANIM_DPI,
+ savefig_kwargs={"facecolor": fig.get_facecolor()})
_extract_poster(webm)
made_gif = False
if make_gif and _LANG == "en":
@@ -14440,6 +14645,350 @@ def update(k: int) -> tuple[Any, ...]:
frames=int(tot_all.shape[1]), gif_fps=8)
+_META_DX = 0.00025
+_META_NY, _META_NX = 4800, 6400
+_META_XL, _META_FACE = 0.45, 0.36
+_META_PITCH, _META_PERIODS = 0.07, 2
+_META_F0 = 2000.0
+
+
+def _metadiffuser_panel_mask(rho: Any) -> None:
+ """Carve the Table-1 metadiffuser (two periods) into a dense slab.
+
+ The slab spans the panel depth L = 2 cm plus a 3 mm back wall under
+ the face line; each well is a vertical slit from the face with its two
+ resonators shelved sideways into the septum, the same layout the
+ to-scale drawing uses (slit at 0.12 d into the cell, lattice a = L/2).
+ """
+ dx, y1 = _META_DX, _META_FACE
+ rows = _METADIFFUSER_T1_ROWS
+ depth, back = 0.02, 0.003
+ rho[round((y1 - depth - back) / dx):round(y1 / dx),
+ round(_META_XL / dx):round((_META_XL + _META_PERIODS * 5
+ * _META_PITCH) / dx)] = 1.2e6
+ lattice = depth / 2
+ for period in range(_META_PERIODS):
+ for n, (h, ln, lc, wn, wc) in enumerate(rows):
+ x0 = _META_XL + (period * 5 + n) * _META_PITCH
+ x_slit = x0 + 0.12 * _META_PITCH
+ c0s, c1s = round(x_slit / dx), round((x_slit + h * 1e-3) / dx)
+ rho[round((y1 - depth) / dx):round(y1 / dx), c0s:c1s] = 1.2
+ for m in range(2):
+ y_m = y1 - depth + (2 - m - 0.5) * lattice
+ x_neck = x_slit + h * 1e-3
+ r0 = round((y_m - 0.5 * wn * 1e-3) / dx)
+ r1 = round((y_m + 0.5 * wn * 1e-3) / dx)
+ rho[r0:r1, c1s:round((x_neck + ln * 1e-3) / dx)] = 1.2
+ r0 = round((y_m - 0.5 * wc * 1e-3) / dx)
+ r1 = round((y_m + 0.5 * wc * 1e-3) / dx)
+ rho[r0:r1, round((x_neck + ln * 1e-3) / dx):
+ round((x_neck + (ln + lc) * 1e-3) / dx)] = 1.2
+
+
+def _meta_qrd_wells() -> list[tuple[float, float, float]]:
+ """QRD well openings (x0, x1, depth) matched to the metadiffuser.
+
+ The same N = 5 residue order as the Table-1 metadiffuser
+ (s = 1, 4, 4, 1, 0), designed at 500 Hz: depths s lambda0 / (2 N) up
+ to 27.4 cm, 65 mm wells split by 5 mm fins on the 7 cm pitch.
+ """
+ unit = (343.0 / 500.0) / 10.0
+ wells: list[tuple[float, float, float]] = []
+ for period in range(_META_PERIODS):
+ for n, s_n in enumerate((1, 4, 4, 1, 0)):
+ x0 = _META_XL + (period * 5 + n) * _META_PITCH
+ wells.append((x0 + 0.005, x0 + _META_PITCH, s_n * unit))
+ return wells
+
+
+def _meta_rho(kind: str) -> Any:
+ """Density map of one run: ``flat``, ``qrd``, ``meta`` or ``ref``."""
+ dx, ny, nx = _META_DX, _META_NY, _META_NX
+ y1 = _META_FACE
+ x_r = _META_XL + _META_PERIODS * 5 * _META_PITCH
+ rho = np.full((ny, nx), 1.2)
+ if kind == "qrd":
+ rho[round(0.05 / dx):round(y1 / dx),
+ round(_META_XL / dx):round(x_r / dx)] = 1.2e6
+ for wx0, wx1, d in _meta_qrd_wells():
+ if d > 0.0:
+ rho[round((y1 - d) / dx):round(y1 / dx),
+ round(wx0 / dx):round(wx1 / dx)] = 1.2
+ elif kind == "meta":
+ _metadiffuser_panel_mask(rho)
+ elif kind == "flat":
+ # A flat rigid slab with the metadiffuser's exact silhouette, so
+ # the fans can only come from the slits, not from the outline.
+ rho[round((y1 - 0.023) / dx):round(y1 / dx),
+ round(_META_XL / dx):round(x_r / dx)] = 1.2e6
+ return rho
+
+
+def _meta_taper() -> Any:
+ """Lateral cosine taper of the incident packet (free-field edges).
+
+ The wavefront is flat over the panels and dies smoothly well before
+ the lateral sponges, so the exterior boundaries only ever absorb
+ outgoing scattered waves and the incident front stays plane: no edge
+ arcs in the total field, and the reference run subtracts exactly.
+ """
+ x = (np.arange(_META_NX) + 0.5) * _META_DX
+ x0, x1, x2, x3 = 0.10, 0.30, 1.30, 1.50
+ w = np.zeros(_META_NX)
+ rise = (x >= x0) & (x < x1)
+ w[rise] = 0.5 - 0.5 * np.cos(np.pi * (x[rise] - x0) / (x1 - x0))
+ w[(x >= x1) & (x <= x2)] = 1.0
+ fall = (x > x2) & (x <= x3)
+ w[fall] = 0.5 + 0.5 * np.cos(np.pi * (x[fall] - x2) / (x3 - x2))
+ return w
+
+
+_META_FRAMES = 600 # 6.13 ms / 10.2 us per frame (49 frames per period)
+
+
+def _meta_pair_worker(kind: str, n_frames: int) -> tuple[Any, Any, Any]:
+ """One panel run in lockstep with its own free-field reference.
+
+ Each worker process pays for its private incident-field run, but the
+ three workers advance in parallel, so the wall time of the cached
+ field computation is about two simulations instead of four.
+ """
+ import fdtd2d
+
+ c0, dx = 343.0, _META_DX
+ y1 = _META_FACE
+ lam0 = c0 / _META_F0
+ # Duration rule: full flight source -> deepest well bottom -> top of
+ # the visible frame, (0.714 + 1.034) m / c * 1.2 = 6.12 ms, sampled at
+ # four times the 12 frames-per-period visual floor (49/period at 2 kHz).
+ every = 33
+ sims = []
+ for rho in (_meta_rho(kind), _meta_rho("ref")):
+ sim = fdtd2d.FDTD2D(c0, dx, rho=rho, shape=(_META_NY, _META_NX),
+ sponge_width=200)
+ sim.add_plane_wave("up", center=0.80, width=0.05, wavelength=lam0)
+ taper = _meta_taper()
+ sim.p *= taper[np.newaxis, :]
+ sim.vy *= taper[np.newaxis, :]
+ sims.append(sim)
+ decay = float(2.0 ** (-sims[0].dt / 0.0015))
+ face_row = round(y1 / dx)
+ trail = np.zeros_like(sims[0].p)
+ tot_frames: list[Any] = []
+ trail_frames: list[Any] = []
+ ts: list[float] = []
+ for _ in range(every * n_frames):
+ for sim in sims:
+ sim.step()
+ scat = sims[0].p - sims[1].p
+ scat[:face_row, :] = 0.0
+ np.maximum(trail * decay, np.abs(scat), out=trail)
+ if sims[0].n % every == 0 and len(ts) < n_frames:
+ tot_frames.append(sims[0].p[::20, ::20].astype(np.float32))
+ trail_frames.append(trail[::20, ::20].astype(np.float32))
+ ts.append(sims[0].time)
+ return np.stack(tot_frames), np.stack(trail_frames), np.asarray(ts)
+
+
+@lru_cache(maxsize=1)
+def _metadiffuser_fields(
+ n_frames: int = _META_FRAMES,
+) -> tuple[Any, Any, Any]:
+ """Plane-wave packet onto a deep QRD vs the 2 cm metadiffuser, cached.
+
+ A 1.6 m x 1.2 m free-field box at dx = 0.25 mm, fine enough to mesh the
+ real millimetre slits, necks and cavities of the Table-1 metadiffuser
+ (density contrast 1e6:1 builds the panels). Each panel (flat control
+ with the metadiffuser's silhouette, deep QRD, metadiffuser) advances
+ in lockstep with a free-field reference inside its own worker process,
+ so ``total - incident`` is exact; a downgoing carrier packet at the
+ 2 kHz evaluation frequency plays the goniometer source. Returns the total
+ frame stacks, the fading scattered-envelope trails [dB] and the frame
+ times. The quantitative far-field comparison lives in the companion
+ polar figure; this clip shows the near fields of the meshed panels.
+ """
+ import multiprocessing as mp
+
+ tot, trail, times = None, None, None
+ try:
+ # The GPU host of .env turns the half-hour CPU computation into a
+ # couple of minutes (the runner falls back by itself, but the CPU
+ # pool below is faster than its single-process local mode).
+ import fdtd_gpu_remote
+
+ fdtd_gpu_remote.load_env()
+ config = fdtd_gpu_remote.RemoteConfig.from_env()
+ use_gpu = bool(config.host) and fdtd_gpu_remote.remote_available(
+ config)
+ except (ImportError, OSError, ValueError):
+ use_gpu = False
+ if use_gpu:
+ import fdtd2d
+
+ # Same duration rule as the CPU worker: 6.12 ms at 49 frames/period.
+ every = 33
+ stride = 20
+ dt = fdtd2d.FDTD2D(343.0, _META_DX, shape=(4, 8)).dt
+ sample_steps = [every * k for k in range(1, n_frames + 1)]
+ lam0 = 343.0 / _META_F0
+ frames: dict[str, Any] = {}
+ for kind in ("flat", "qrd", "meta", "ref"):
+ job = fdtd_gpu_remote.build_job(
+ 343.0, _META_DX, steps=every * n_frames,
+ sample_steps=sample_steps, shape=(_META_NY, _META_NX),
+ rho=_meta_rho(kind), sponge_width=200,
+ plane_waves=[{"direction": "up", "center": 0.80,
+ "width": 0.05, "wavelength": lam0}],
+ init_scale_x=_meta_taper(),
+ sample_stride=stride, sample_dtype="float32",
+ )
+ # The four 19 800-step field jobs run ~10-12 min each on the
+ # GPU host; give each submit an hour before falling back.
+ frames[kind] = fdtd_gpu_remote.submit(
+ job, config, timeout=3600.0)["frames"]
+ face_row = round(_META_FACE / _META_DX) // stride
+ decay = float(2.0 ** (-(every * dt) / 0.0015))
+ tot_list, trail_list = [], []
+ for kind in ("flat", "qrd", "meta"):
+ scat = frames[kind] - frames["ref"]
+ scat[:, :face_row, :] = 0.0
+ running = np.zeros_like(scat[0])
+ history = np.empty_like(scat)
+ for k in range(scat.shape[0]):
+ np.maximum(running * decay, np.abs(scat[k]), out=running)
+ history[k] = running
+ tot_list.append(frames[kind].astype(np.float32))
+ trail_list.append(history)
+ tot = np.stack(tot_list)
+ trail = np.stack(trail_list)
+ times = np.asarray(sample_steps, dtype=np.float64) * dt
+ else:
+ ctx = mp.get_context("spawn")
+ with ctx.Pool(processes=3) as pool:
+ parts = pool.starmap(
+ _meta_pair_worker,
+ [(kind, n_frames) for kind in ("flat", "qrd", "meta")],
+ )
+ tot = np.stack([part[0] for part in parts])
+ trail = np.stack([part[1] for part in parts])
+ times = parts[0][2]
+ ref = float(trail[:, trail.shape[1] // 3:].max()) or 1.0
+ with np.errstate(divide="ignore"):
+ trail_db = 20.0 * np.log10(trail / ref)
+ trail_db = np.clip(trail_db, -30.0, 0.0).astype(np.float32)
+ return tot, trail_db, times
+
+
+def animate_fdtd_metadiffuser(output_dir: str) -> None:
+ """The 27 cm deep Schroeder QRD vs the 2 cm metadiffuser that mimics
+ it, next to a flat control slab (2D FDTD at 0.25 mm, real slits, necks
+ and cavities meshed): the same 2 kHz wavefront leaves the same kind of
+ scattered fan, from a panel 13.7 times thinner."""
+ from matplotlib import patheffects
+ from matplotlib.patches import Polygon, Rectangle
+
+ T = _translate_str
+ outline = [patheffects.withStroke(linewidth=2.0, foreground="white")]
+ tot_all, trail_db, times = _metadiffuser_fields()
+ y1 = _META_FACE
+ x_l = _META_XL
+ x_r = x_l + _META_PERIODS * 5 * _META_PITCH
+ vmax = float(np.quantile(np.abs(tot_all[:, 0]), 0.999))
+
+ fig = _anim_figure()
+ fig.suptitle(T("Schroeder diffuser vs metadiffuser (2D FDTD)"),
+ fontweight="bold")
+ gs = fig.add_gridspec(2, 3)
+ titles = [T("Flat rigid panel"), T("QRD, wells down to 27 cm"),
+ T("Metadiffuser, 2 cm panel")]
+ qrd_poly: list[tuple[float, float]] = [(x_l, 0.05), (x_l, y1)]
+ for wx0, wx1, d in _meta_qrd_wells():
+ qrd_poly += [(wx0, y1), (wx0, y1 - d), (wx1, y1 - d), (wx1, y1)]
+ qrd_poly += [(x_r, y1), (x_r, 0.05)]
+ xc = 0.5 * (x_l + x_r)
+
+ ims: list[Any] = []
+ d_txts: list[Any] = []
+ for col in range(3):
+ ax_t = fig.add_subplot(gs[0, col])
+ ax_s = fig.add_subplot(gs[1, col])
+ im_t = ax_t.imshow(tot_all[col][0], origin="lower",
+ extent=(0.0, 1.6, 0.0, 1.2), cmap="RdBu_r",
+ vmin=-vmax, vmax=vmax, interpolation="bilinear")
+ im_s = ax_s.imshow(trail_db[col][0], origin="lower",
+ extent=(0.0, 1.6, 0.0, 1.2), cmap="magma",
+ vmin=-30.0, vmax=0.0, interpolation="bilinear")
+ ax_t.set_title(titles[col], fontsize=10, fontweight="bold")
+ for ax in (ax_t, ax_s):
+ ax.grid(False)
+ if col == 1:
+ ax.add_patch(Polygon(qrd_poly, closed=True,
+ facecolor=COLOR_GRID,
+ edgecolor=COLOR_FG, lw=0.8))
+ else:
+ ax.add_patch(Rectangle((x_l, y1 - 0.023), x_r - x_l, 0.023,
+ facecolor=COLOR_GRID,
+ edgecolor=COLOR_FG, lw=0.8))
+ ax.set_xlim(0.06, 1.54)
+ ax.set_ylim(0.0, 1.12)
+ ax.tick_params(labelsize=7)
+ ax_t.tick_params(labelbottom=False)
+ ax_s.set_xlabel("x [m]", fontsize=8)
+ if col == 1:
+ ax_t.text(xc, 0.97, T("incident plane wavefront"), ha="center",
+ va="bottom", color="black", fontsize=7.5,
+ path_effects=outline)
+ ax_t.annotate("", xy=(xc, 0.83), xytext=(xc, 0.955),
+ arrowprops={"arrowstyle": "-|>", "color": "black",
+ "lw": 1.2})
+ d_txt = ax_s.text(xc, 1.03, "", ha="center", va="top",
+ color="white", fontsize=7.5, fontweight="bold")
+ if col == 0:
+ ax_t.set_ylabel(T("sound field p"), fontsize=9)
+ ax_s.set_ylabel(T("scattered field (total − incident)"),
+ fontsize=8)
+ else:
+ ax_t.tick_params(labelleft=False)
+ ax_s.tick_params(labelleft=False)
+ if col == 1:
+ ax_t.annotate("", xy=(x_r + 0.045, y1 - 0.274),
+ xytext=(x_r + 0.045, y1),
+ arrowprops={"arrowstyle": "-", "color": COLOR_FG,
+ "lw": 1.6})
+ ax_t.text(x_r + 0.07, y1 - 0.14, "27 cm", ha="left",
+ va="center", fontsize=7, color=COLOR_FG)
+ if col == 2:
+ ax_t.text(xc, 0.06, T("real slits and resonators meshed at "
+ "0.25 mm"), ha="center", va="bottom",
+ color="black", fontsize=6.5, path_effects=outline)
+ ax_t.annotate("", xy=(x_r + 0.045, y1 - 0.023),
+ xytext=(x_r + 0.045, y1),
+ arrowprops={"arrowstyle": "-", "color": COLOR_FG,
+ "lw": 1.6})
+ ax_t.text(x_r + 0.07, y1 - 0.012, "2 cm", ha="left",
+ va="center", fontsize=7, color=COLOR_FG)
+ ims += [im_t, im_s]
+ d_txts.append(d_txt)
+ t_txt = fig.text(0.985, 0.02, "", ha="right", va="bottom",
+ family="monospace", fontsize=10, color=COLOR_FG)
+ reveal = int(0.8 * tot_all.shape[1])
+
+ verdicts = [T("a collimated specular beam"),
+ T("a wide scattered fan"),
+ T("the same fan, from 2 cm")]
+
+ def update(k: int) -> tuple[Any, ...]:
+ for col in range(3):
+ ims[2 * col].set_data(tot_all[col][k])
+ ims[2 * col + 1].set_data(trail_db[col][k])
+ d_txts[col].set_text(verdicts[col] if k >= reveal else "")
+ t_txt.set_text(T(f"t = {times[k] * 1000.0:4.2f} ms"))
+ return (*ims, *d_txts, t_txt)
+
+ _render_clip(fig, update, output_dir, "anim_fdtd_metadiffuser",
+ frames=int(tot_all.shape[1]), gif_fps=8)
+
+
def animate_standing_wave_tube(output_dir: str) -> None:
"""ISO 10534-2 impedance tube: the incident and reflected waves travel
inside a drawn tube and their sum forms the standing-wave envelope; a
@@ -15493,6 +16042,7 @@ def update(kf: int) -> tuple[Any, ...]:
"anim_fdtd_ground_effect": animate_fdtd_ground_effect,
"anim_fdtd_ducting": animate_fdtd_ducting,
"anim_fdtd_diffusion": animate_fdtd_diffusion,
+ "anim_fdtd_metadiffuser": animate_fdtd_metadiffuser,
"anim_fdtd_impedance_tube": animate_fdtd_impedance_tube,
"anim_fdtd_transmission_tube": animate_fdtd_transmission_tube,
"anim_standing_wave_tube": animate_standing_wave_tube,
@@ -15707,6 +16257,7 @@ def _generate_figures_parallel(
"anim_fdtd_ground_effect": 300.0,
"anim_fdtd_ducting": 280.0,
"anim_fdtd_diffusion": 260.0,
+ "anim_fdtd_metadiffuser": 540.0,
"anim_standing_wave_tube": 130.0,
"anim_sweep_deconvolution": 90.0,
"anim_power_two_rooms": 80.0,
diff --git a/site/src/content/docs/es/guides/surface-scattering.mdx b/site/src/content/docs/es/guides/surface-scattering.mdx
index 056edd603..bdffc94fd 100644
--- a/site/src/content/docs/es/guides/surface-scattering.mdx
+++ b/site/src/content/docs/es/guides/surface-scattering.mdx
@@ -10,6 +10,15 @@ references:
publisher: "CRC Press"
doi: "10.1201/9781315369211"
note: "ISBN 978-1-4987-4099-9. La monografía de referencia sobre teoría y diseño de difusores, de los autores del método del coeficiente de difusión de ISO 17497-2: la distinción dispersión-difusión, los montajes de medida y la guía de diseño que esta página condensa. El Apéndice B (tabla de coeficientes de difusión normalizados, pp. 481-485) es el anclaje BEM publicado del modelo de predicción de difusores de esta página."
+ - type: article
+ authors: ["Jiménez, N.", "Cox, T. J.", "Romero-García, V.", "Groby, J.-P."]
+ year: 2017
+ title: "Metadiffusers: Deep-subwavelength sound diffusers"
+ container: "Scientific Reports"
+ volume: "7"
+ pages: "5389"
+ doi: "10.1038/s41598-017-05710-5"
+ note: "El modelo de metadifusor implementado aquí: rendijas cargadas con resonadores de Helmholtz reproducen perfiles de fase de Schroeder y secuencias ternarias con paneles de 1/46 a 1/20 de la longitud de onda de diseño."
- type: article
authors: ["Hargreaves, T. J.", "Cox, T. J.", "Lam, Y. W.", "D'Antonio, P."]
year: 2000
@@ -586,6 +595,134 @@ plt.show()
+### Metadifusores: difusores de Schroeder en sublongitud de onda profunda
+
+Un metadifusor sustituye los pozos profundos de un difusor de Schroeder por
+rendijas finas cargadas con resonadores de Helmholtz (Jiménez, Cox,
+Romero-García y Groby, 2017). Por debajo de su resonancia, los resonadores
+ralentizan el sonido dentro de cada rendija, de modo que un panel de pocos
+centímetros alcanza las fases de reflexión que a una red de fase clásica le
+exigen decenas de centímetros de profundidad; y llevar una rendija al
+acoplamiento crítico añade un estado perfectamente absorbente, el `0` que
+requieren las secuencias ternarias. `metadiffuser_reflection` ejecuta la
+cadena de matrices de transferencia del
+[absorbedor de sonido lento](/phonometry/es/guides/materials/) una vez por
+pozo (resonadores bidimensionales, pérdidas viscotérmicas y correcciones de
+extremo incluidas) y devuelve la reflexión compleja por pozo $R_n(f)$;
+`metadiffuser_polar_response` y `metadiffuser_diffusion_spectrum` reducen
+ese perfil espacial con el mismo campo lejano de Fraunhofer y el mismo
+coeficiente ISO 17497-2 que los diseños clásicos de arriba.
+
+
+Ver el código de esta figura
+
+```python
+import numpy as np
+from phonometry import (
+ HelmholtzResonator,
+ MetadiffuserWell,
+ materials,
+ metadiffuser_polar_response,
+)
+
+mm = 1e-3
+filas = [
+ (14.7, 13.0, 16.4, 6.2, 9.0),
+ (30.9, 9.1, 4.3, 3.5, 9.0),
+ (30.9, 9.1, 4.3, 3.5, 9.0),
+ (15.7, 13.3, 17.0, 6.3, 9.0),
+ (20.3, 18.0, 20.7, 3.2, 9.0),
+]
+pozos = [
+ MetadiffuserWell(
+ h * mm,
+ 2 * (HelmholtzResonator(ln * mm, wn * mm, lc * mm, wc * mm),),
+ )
+ for h, ln, lc, wn, wc in filas
+]
+
+# El panel metadifusor y el QRD al que se ajustó, a 2 kHz y con seis
+# repeticiones del periodo.
+meta = metadiffuser_polar_response(2000.0, pozos, depth=0.02, period=0.07,
+ periods=6)
+secuencia = np.roll(materials.quadratic_residue_sequence(5), -1)
+profundidades = secuencia * (343.0 / 500.0) / (2 * 5)
+qrd = materials.predict_diffuser_polar_response(
+ 0.07, 2000.0, depths=profundidades, periods=6, include_obliquity=False,
+)
+
+ax = meta.plot(marker="", linewidth=2.2, language="es",
+ label="Metadifusor, panel de 2 cm")
+qrd.plot(ax=ax, marker="", linewidth=1.6, linestyle="--", language="es",
+ label="QRD, pozos de hasta 27,4 cm")
+ax.legend(loc="lower center")
+```
+
+
+
+La superposición de campo lejano de arriba es el resumen en frecuencia; la
+animación FDTD de abajo malla ambos paneles de verdad (el metadifusor a
+0,25 mm, con sus rendijas, cuellos y cavidades) y lanza el mismo frente de
+2 kHz contra ellos: el QRD de 27 cm y el panel de 2 cm devuelven abanicos
+dispersados casi idénticos.
+
+
+Show the code for this figure
+
+```python
+import numpy as np
+from phonometry import (
+ HelmholtzResonator,
+ MetadiffuserWell,
+ materials,
+ metadiffuser_polar_response,
+)
+
+mm = 1e-3
+rows = [
+ (14.7, 13.0, 16.4, 6.2, 9.0),
+ (30.9, 9.1, 4.3, 3.5, 9.0),
+ (30.9, 9.1, 4.3, 3.5, 9.0),
+ (15.7, 13.3, 17.0, 6.3, 9.0),
+ (20.3, 18.0, 20.7, 3.2, 9.0),
+]
+wells = [
+ MetadiffuserWell(
+ h * mm,
+ 2 * (HelmholtzResonator(ln * mm, wn * mm, lc * mm, wc * mm),),
+ )
+ for h, ln, lc, wn, wc in rows
+]
+
+# The metadiffuser panel and the QRD it was tuned to at 2 kHz, both with
+# six repetitions of the period.
+meta = metadiffuser_polar_response(2000.0, wells, depth=0.02, period=0.07,
+ periods=6)
+sequence = np.roll(materials.quadratic_residue_sequence(5), -1)
+depths = sequence * (343.0 / 500.0) / (2 * 5)
+qrd = materials.predict_diffuser_polar_response(
+ 0.07, 2000.0, depths=depths, periods=6, include_obliquity=False,
+)
+
+ax = meta.plot(marker="", linewidth=2.2, label="Metadiffuser, panel 2 cm")
+qrd.plot(ax=ax, marker="", linewidth=1.6, linestyle="--",
+ label="QRD, wells up to 27.4 cm")
+ax.legend(loc="lower center")
+```
+
+
+
+The far-field overlay above is the frequency-domain summary; the FDTD
+animation below meshes both panels for real (the metadiffuser at 0.25 mm, slits,
+necks and cavities included) and lets the same 2 kHz wavefront hit them: the
+27 cm QRD and the 2 cm panel throw out nearly the same scattered fan.
+
+ 
