Lo mejor de dos mundos: DGTD

FDTD es intuitivo y broadband, pero su rejilla cartesiana escalonea las curvas. FEM tiene mallas elegantes que siguen la geometría, pero trabaja en dominio frecuencial — una simulación por frecuencia. ¿Y si pudieras tener ambas cosas? Mallas triangulares adaptativas Y evolución temporal. Eso es DGTD.

El problema del escalonado

En el artículo 04 vimos FDTD corriendo en tiempo real. Funciona, pero tiene un problema fundamental: la rejilla de Yee es cartesiana. Cuando el objeto tiene curvas — un cilindro, una esfera, un resonador de anillo — la rejilla lo aproxima con escalones. Eso introduce errores de fase que limitan la precisión a segundo orden, por mucho que refines la rejilla.

Para un resonador óptico de alto $Q$ , esos errores de fase son fatales. La frecuencia de resonancia converge lentamente, y el factor de calidad $Q$ puede estar mal por órdenes de magnitud incluso con rejillas finas (Niegemann et al., J. Opt. A, 2009).

Compara las dos aproximaciones del mismo cilindro:

Explorar

Resolución 10

Mismo cilindro, misma resolución. FDTD escalonea las curvas (rojo). DGTD las sigue fielmente (verde).

A la izquierda, FDTD: la línea roja punteada es el contorno real. Los bloques azules son las celdas que la rejilla clasifica como «dentro del cilindro». Las escaleras son inevitables. A la derecha, DGTD: los triángulos se adaptan al contorno y lo siguen sin escalones. Aumenta la resolución y verás cómo el FDTD sigue escalonando mientras el DGTD refina suavemente.

¿Cómo funciona?

DGTD toma prestado de tres tradiciones:

De FEM: mallas no estructuradas (triángulos en 2D, tetraedros en 3D) con funciones base polinómicas de orden $p$ en cada elemento.
De FDTD: evolución temporal explícita (típicamente Runge-Kutta de orden 4). No hay sistemas lineales grandes que resolver en cada paso de tiempo.
De volúmenes finitos: el acoplamiento entre elementos vecinos se maneja mediante flujos numéricos en las interfaces. Como las funciones base son locales a cada elemento (de ahí «discontinuas»), los campos pueden ser diferentes a cada lado de una interfaz — la información se intercambia solo a través del flujo.

Matemáticamente, en cada elemento $\Delta$ se proyectan las ecuaciones de Maxwell sobre funciones test locales $L_i$ :

\int_\Delta \left(\mathscr{Q} \cdot \partial_t \mathbf{q}_N + \nabla \cdot \mathbf{F}(\mathbf{q}_N)\right) \cdot L_i \, d^3r = 0

donde $\mathbf{q} = (E, H)^T$ es el vector de estado y $\mathbf{F}$ es el flujo de Maxwell. La clave: esta ecuación es completamente local — solo involucra campos dentro del elemento. El acoplamiento con los vecinos entra como un término de flujo en la frontera del elemento, que reemplaza $\nabla \cdot \mathbf{F}$ por un flujo numérico (típicamente upwind) evaluado en las caras compartidas.

¿Por qué «discontinuo»?

En FEM convencional, las funciones base son continuas entre elementos — los campos se comparten en los nodos. Esto crea un sistema de ecuaciones global (matriz dispersa grande) que hay que resolver.

En DG, las funciones base son locales a cada elemento. Los campos pueden saltar en la frontera entre elementos — de ahí «discontinuo». Esto hace que la matriz de masas sea bloque-diagonal (un bloque por elemento), trivialmente invertible. El resultado: stepping temporal explícito sin resolver sistemas lineales grandes.

El precio: necesitas un mecanismo de acoplamiento (el flujo numérico) para que la información pase de un elemento a otro. Pero ese flujo es local — solo involucra los dos elementos que comparten la cara.

Convergencia exponencial

La ventaja más espectacular de DGTD sobre FDTD es la convergencia. FDTD converge como $O(h^2)$ — segundo orden. Duplicar la resolución (y multiplicar el coste por 16 en 3D) solo gana un factor 4 en precisión.

DGTD con polinomios de orden $p$ converge como $O(h^{p+1})$ — y para problemas suaves, la convergencia es exponencial en $p$ . Esto significa que subir el orden polinómico de 3 a 4 puede mejorar la precisión más que duplicar la malla.

El paper de Niegemann et al. (2009) lo demuestra con un resonador cilíndrico: FDTD necesita ~128 celdas/radio para bajar el error en $Q$ al 1%. DGTD con orden $p = 3$ lo consigue con ~8 elementos/radio — usando una fracción de la memoria.

CFL en mallas no estructuradas

DGTD, como FDTD, es un método explícito — y tiene su propia condición CFL. Pero en mallas no estructuradas, la condición es más restrictiva:

\Delta t \leq C \, \frac{h_{\min}}{c \, p^2}

donde $h_{\min}$ es el tamaño del elemento más pequeño, $p$ el orden polinómico, y $C$ una constante que depende de la geometría del elemento. El factor $p^2$ es el precio de la convergencia de alto orden: a mayor $p$ , paso temporal más pequeño.

El problema: si la malla tiene unos pocos elementos muy pequeños (cerca de una punta o una rendija), esos elementos dictan el paso temporal de toda la simulación. La solución es el local time stepping (LTS): cada elemento avanza con su propio $\Delta t$ , proporcional a su tamaño. Los elementos grandes dan pasos más largos; los pequeños, más cortos. Esto puede acelerar la simulación 10–100× en problemas con mallas muy heterogéneas.

¿Cuándo usar DGTD?

Resonadores de alto Q: microdiscos, anillos, cavidades fotónicas. FDTD no puede capturar Q > 10⁴ de forma fiable; DGTD sí.
Geometrías curvas: nanopartículas esféricas o elipsoidales, guías de onda curvas, lentes. DGTD las discretiza sin escaleras.
Problemas multiscala: la malla puede ser gruesa en el vacío y ultra-fina cerca de puntas, rendijas o interfaces. La adaptatividad h-p (refinar malla Y subir orden) es natural.
Broadband: como FDTD, un pulso corto te da todo el espectro. Pero con mejor precisión.

El precio: DGTD es más complejo de implementar que FDTD (generación de malla, funciones base de alto orden, flujos numéricos). Pero las implementaciones open-source (la librería de Hesthaven & Warburton, Nodal DG; el código del grupo de Busch, BEAST) y los solvers comerciales lo han democratizado.

DGTD en la práctica

El grupo de Kurt Busch (Humboldt-Universität zu Berlin) ha sido pionero en aplicar DGTD a nanofotónica. Sus trabajos demuestran ventajas claras en:

Resonadores de cristal fotónico: cavidades con $Q > 10^6$ , donde FDTD no puede converger (Busch et al., Laser & Photonics Reviews, 2011).
Plasmónica: nanotriángulos, bow-ties, nanoantenas con puntas afiladas. La malla refina en la punta sin penalizar el paso temporal global (con LTS).
Metamateriales: celdas unitarias con geometrías complejas (split-ring resonators, fish-nets) donde las escalas van de 10 nm a varios μm.

El paisaje completo

Con DGTD cerramos un panorama que cubre todas las necesidades de la simulación en fotónica:

TMM: 1D, instantáneo — capas planas.
Mie/Rayleigh: analítico — esferas y partículas pequeñas.
FDTD: tiempo, rejilla cartesiana — rápido, intuitivo.
DDA: dipolos en espacio abierto — nanopartículas arbitrarias.
FEM: frecuencia, mallas adaptativas — cavidades y guías.
BEM: superficie, espacio abierto — plasmónica, campos cercanos.
DGTD: tiempo + mallas adaptativas — resonadores de alto Q, multiscala.

No hay un método universal. Pero con estos siete (+1), tienes las herramientas para atacar prácticamente cualquier problema de nanofotónica. El arte — como siempre — está en elegir el correcto.

Ejercicios

Ejercicio 1

Un resonador de microdisco tiene un factor de calidad $Q = 10^5$ . Eso significa que la luz da $Q/2\pi \approx 16\,000$ vueltas dentro del disco antes de disiparse. Si el disco tiene radio $R = 5$ μm y opera a $\lambda = 1.55$ μm, ¿cuántos periodos ópticos dura la simulación temporal? ¿Y cuántos pasos de tiempo si $\Delta t \approx T/20$ (donde $T$ es el periodo óptico)?

Solución

Periodo óptico: $T = \lambda/c = 1.55 \times 10^{-6} / 3 \times 10^8 \approx 5.2$ fs. Tiempo de vida del fotón en el resonador: $\tau = Q \cdot T / 2\pi = 10^5 \times 5.2 / 6.28 \approx 83\,000$ fs = 83 ps.

Periodos ópticos: $\tau/T \approx 16\,000$ . Pasos de tiempo: $16\,000 \times 20 = 320\,000$ . Con FDTD (rejilla fina para reducir staircasing) cada paso es caro. Con DGTD (malla gruesa + alto orden), cada paso cuesta menos y la convergencia es mejor. Pero 320.000 pasos son 320.000 pasos — por eso los resonadores de alto Q son un desafío computacional incluso con DGTD.

Ejercicio 2

Compara los costes de convergencia. FDTD tiene error $O(h^2)$ ; DGTD con orden $p$ tiene error $O(h^{p+1})$ . Si quieres reducir el error por un factor 100 (dos órdenes de magnitud):

Con FDTD, ¿por cuánto necesitas reducir $h$ ? ¿Cuánto sube el coste en 3D?
Con DGTD a $p = 3$ , ¿por cuánto necesitas reducir $h$ ? ¿Cuánto sube el coste?

Solución

FDTD: error $\propto h^2$ . Reducir 100×: $h \to h/10$ . Coste 3D: $(10)^3 \times 10 = 10^4\times$ más ( $N^3$ celdas × $N$ pasos por CFL).

DGTD p=3: error $\propto h^4$ . Reducir 100×: $h \to h/\sqrt[4]{100} \approx h/3.16$ . Coste 3D: $(3.16)^3 \times 3.16 \approx 100\times$ más.

10.000× vs 100× — dos órdenes de magnitud de diferencia. Esto explica por qué DGTD es indispensable para problemas que requieren alta precisión: resonadores de alto Q, propagación larga, o geometrías con features sub-λ.