Derivación del método FDTD.
1. Introducción
Esta página presenta la derivación matemática detallada de las principales ecuaciones de actualización del método de diferencias finitas en el dominio del tiempo (FDTD) que se utilizan para resolver las ecuaciones de Maxwell en el dominio del tiempo. El propósito de esta página es diferente al de la página principal sobre la teoría. El objetivo aquí no es ofrecer una descripción general simplificada, sino mostrar explícitamente cómo las ecuaciones continuas del rotacional se transforman en las expresiones algebraicas discretas que se utilizan en un programa informático. Esta es también la parte del método FDTD en la que es más fácil cometer errores en el código: errores de signo, índices mal colocados o un escalonamiento temporal inconsistente producirán inmediatamente resultados inestables o poco realistas.
Por lo tanto, partimos de la ley de Ampère y la ley de Faraday en forma diferencial, las desarrollamos en componentes cartesianos y, a continuación, las adaptamos a un problema óptico bidimensional. Las derivadas espaciales se sustituyen por diferencias finitas centradas, las derivadas temporales se discretizan utilizando el esquema estándar de «salto de rana» y las expresiones resultantes se reorganizan para obtener fórmulas explícitas de actualización para los componentes del campo eléctrico y magnético.
Para que la derivación sea más fácil de entender, reducimos aquí las ecuaciones a un caso bidimensional en el que \(\partial/\partial y = 0\) y la conductividad magnética se establecen en cero, \(\sigma_m = 0\). Esto simplifica considerablemente el álgebra y también muestra explícitamente cómo se lleva a cabo la reducción de dimensiones en la práctica. Sin embargo, en el solucionador OghmaNano FDTD completo, las ecuaciones electromagnéticas se resuelven en tres dimensiones. Por lo tanto, la derivación bidimensional que se presenta aquí tiene como objetivo ser un ejemplo claro y práctico del método, en lugar de una afirmación de que el solucionador subyacente está limitado a 2D.
A lo largo de esta página, se evalúan los componentes del campo eléctrico en intervalos de tiempo enteros \(n\Delta t\), mientras que los componentes del campo magnético se evalúan en intervalos de tiempo de media unidad \(\left(n+\tfrac{1}{2}\right)\Delta t\). Este es el esquema estándar de la malla de Yee utilizado en el método FDTD. A diferencia de la página de resumen más breve, aquí se indica de forma más explícita la indexación para que la derivación pueda leerse tanto como una guía de implementación como una nota matemática.
2. Ecuaciones iniciales
Partimos de las ecuaciones de Maxwell en forma vectorial para un medio lineal con conductividad eléctrica \(\sigma\), permitividad eléctrica \(\epsilon\), permeabilidad magnética \(\mu\) y, para mayor generalidad, conductividad magnética \(\sigma_m\):
\[ \nabla \times \boldsymbol{H} = \sigma \boldsymbol{E} + \epsilon \frac{\partial \boldsymbol{E}}{\partial t}, \qquad \nabla \times \boldsymbol{E} = - \sigma_m \boldsymbol{H} - \mu \frac{\partial \boldsymbol{H}}{\partial t}. \]
En forma de matriz, las curvas son
\[ \nabla \times \boldsymbol{H} = \begin{vmatrix} \hat{\boldsymbol{x}} & \hat{\boldsymbol{y}} & \hat{\boldsymbol{z}} \\ \dfrac{\partial}{\partial x} & \dfrac{\partial}{\partial y} & \dfrac{\partial}{\partial z} \\ H_x & H_y & H_z \end{vmatrix}, \qquad \nabla \times \boldsymbol{E} = \begin{vmatrix} \hat{\boldsymbol{x}} & \hat{\boldsymbol{y}} & \hat{\boldsymbol{z}} \\ \dfrac{\partial}{\partial x} & \dfrac{\partial}{\partial y} & \dfrac{\partial}{\partial z} \\ E_x & E_y & E_z \end{vmatrix}. \]
Estas son las ecuaciones en tiempo continuo y espacio continuo. El objetivo del resto de la página es transformarlas en un esquema explícito de diferencias finitas que pueda utilizarse para avanzar en el tiempo.
3. Forma vectorial de las ecuaciones del rotacional.
Al expandir la ley de Ampère en componentes cartesianos, obtenemos:
\[ \sigma E_x + \epsilon \frac{\partial E_x}{\partial t} = \frac{\partial H_z}{\partial y} - \frac{\partial H_y}{\partial z}, \]
\[ \sigma E_y + \epsilon \frac{\partial E_y}{\partial t} = \frac{\partial H_x}{\partial z} - \frac{\partial H_z}{\partial x}, \]
\[ \sigma E_z + \epsilon \frac{\partial E_z}{\partial t} = \frac{\partial H_y}{\partial x} - \frac{\partial H_x}{\partial y}. \]
Al ampliar la ley de Faraday, obtenemos:
\[ -\sigma_m H_x - \mu \frac{\partial H_x}{\partial t} = \frac{\partial E_z}{\partial y} - \frac{\partial E_y}{\partial z}, \]
\[ -\sigma_m H_y - \mu \frac{\partial H_y}{\partial t} = \frac{\partial E_x}{\partial z} - \frac{\partial E_z}{\partial x}, \]
\[ -\sigma_m H_z - \mu \frac{\partial H_z}{\partial t} = \frac{\partial E_y}{\partial x} - \frac{\partial E_x}{\partial y}. \]
En esta etapa, todavía no se ha supuesto nada sobre la dimensionalidad. Las expresiones anteriores son simplemente la forma en componentes de las ecuaciones completas del rotacional tridimensional.
4. Reducción bidimensional
Para facilitar la derivación, reducimos ahora el problema a dos dimensiones, asumiendo que no hay variación en la dirección \(y\), por lo que \(\partial / \partial y = 0\). Este paso se realiza únicamente para simplificar el álgebra y facilitar la comprensión de la estructura de las ecuaciones de actualización. El mismo procedimiento puede extenderse directamente a tres dimensiones y, en la práctica, OghmaNano resuelve la forma completa en 3D de las ecuaciones de Maxwell.
Partiendo de esta premisa, las ecuaciones del campo eléctrico se simplifican a:
\[ \sigma E_x + \epsilon \frac{\partial E_x}{\partial t} = - \frac{\partial H_y}{\partial z}, \]
\[ \sigma E_y + \epsilon \frac{\partial E_y}{\partial t} = - \frac{\partial H_z}{\partial x} + \frac{\partial H_x}{\partial z}, \]
\[ \sigma E_z + \epsilon \frac{\partial E_z}{\partial t} = \frac{\partial H_y}{\partial x}. \]
Si también fijamos la conductividad magnética en cero, \(\sigma_m = 0\), las ecuaciones del campo magnético se convierten en:
\[ \mu \frac{\partial H_x}{\partial t} = \frac{\partial E_y}{\partial z}, \]
\[ \mu \frac{\partial H_y}{\partial t} = \frac{\partial E_z}{\partial x} - \frac{\partial E_x}{\partial z}, \]
\[ \mu \frac{\partial H_z}{\partial t} = - \frac{\partial E_y}{\partial x}. \]
Estas son las ecuaciones continuas que ahora se discretizarán. Es importante señalar explícitamente que se trata de una reducción bidimensional del problema tridimensional completo. El objetivo no es afirmar que todos los métodos FDTD son intrínsecamente bidimensionales, sino obtener ecuaciones de actualización claras para el caso común en el que se suprime una dirección.
5. Cuadrícula de Yee, indexación y notación.
Para resolver numéricamente las ecuaciones anteriores, tanto el espacio como el tiempo se discretizan. Supongamos que los intervalos de la malla computacional en las direcciones \(x\) y \(z\) son \(\Delta x\) y \(\Delta z\), respectivamente, y que el intervalo de tiempo es \(\Delta t\). Utilizamos los índices enteros \(i\) y \(k\) para las direcciones \(x\) y \(z\):
\[ x_i = i \Delta x, \qquad z_k = k \Delta z. \]
El tiempo se mide mediante el número entero \(n\), de modo que
\[ t^n = n \Delta t. \]
En el esquema de Yee, los campos eléctricos y magnéticos no se almacenan en las mismas ubicaciones, ni en el espacio ni en el tiempo. En cambio, se alternan. Los componentes del campo eléctrico se evalúan en niveles de tiempo enteros \(n\Delta t\), mientras que los componentes del campo magnético se evalúan en niveles de tiempo semienteros \(\left(n+\tfrac{1}{2}\right)\Delta t\). Esta es la base de la estructura de «salto» de FDTD.
Desde el punto de vista espacial, los componentes del campo también están desplazados de manera que cada rotacional discreto centrado esté definido de forma natural entre las ubicaciones de la cuadrícula adyacentes. En una implementación completa del código, cada componente del campo lleva su propio par de índices espaciales, por ejemplo, \(E_x^n(i,k)\), \(H_y^{n+\frac{1}{2}}(i+\tfrac{1}{2},k+\tfrac{1}{2})\), y así sucesivamente. Para que las expresiones algebraicas sean más fáciles de leer, algunos de los índices repetidos, que son meros «observadores», no siempre se escriben explícitamente en los pasos siguientes. Sin embargo, la interpretación siempre es la misma:
- \(i\) indexa la dirección \(x\).
- \(k\) indexa la dirección \(z\).
- \(n\) indica los niveles temporales del campo eléctrico, expresados en números enteros, y
- \(n+\tfrac{1}{2}\) indica los niveles temporales del campo magnético en intervalos de medio tono.
Por ejemplo, la aproximación de diferencias finitas centrada para \(\partial H_y / \partial z\) es:
\[ \frac{\partial H_y}{\partial z} \approx \frac{ H_y^{n+\frac{1}{2}}\!\left(i,k+\frac{1}{2}\right) - H_y^{n+\frac{1}{2}}\!\left(i,k-\frac{1}{2}\right) }{\Delta z}, \]
y, de manera similar, la aproximación centrada a \(\partial H_z / \partial x\) es
\[ \frac{\partial H_z}{\partial x} \approx \frac{ H_z^{n+\frac{1}{2}}\!\left(i+\frac{1}{2},k\right) - H_z^{n+\frac{1}{2}}\!\left(i-\frac{1}{2},k\right) }{\Delta x}. \]
Una vez que se explicita esta notación, la derivación que se presenta a continuación puede interpretarse sin ambigüedades como un esquema adecuado para actualizar la malla de Yee.
6. Obtención de las ecuaciones para la actualización del campo eléctrico.
El campo eléctrico se calcula en intervalos de tiempo enteros \(n\Delta t\), mientras que el campo magnético se calcula en intervalos de tiempo de media unidad \(\left(n+\tfrac{1}{2}\right)\Delta t\). Para obtener las ecuaciones de actualización, la derivada temporal se aproxima mediante una diferencia finita centrada, y el término conductor se calcula utilizando el promedio de los campos eléctricos antiguos y nuevos. Esto proporciona un tratamiento numéricamente estable en forma de trapecio para la conductividad.
Este paso de promediación es importante. Si el término de conductividad se tratara de forma inconsistente, los coeficientes de actualización ya no tendrían la forma estándar atenuada y se perdería la correspondencia con las ecuaciones FDTD habituales.
6.1 Actualizar la ecuación para \(E_x\).
Comienza desde
\[ \sigma E_x + \epsilon \frac{\partial E_x}{\partial t} = - \frac{\partial H_y}{\partial z}. \]
Evalúe el término de conductividad utilizando el promedio de \(E_x^{n+1}\) y \(E_x^n\), y evalúe la derivada temporal utilizando una diferencia finita centrada:
\[ \sigma \frac{E_x^{n+1}(i,k) + E_x^n(i,k)}{2} + \epsilon \frac{E_x^{n+1}(i,k) - E_x^n(i,k)}{\Delta t} = - \frac{ H_y^{n+\frac{1}{2}}\!\left(i,k+\frac{1}{2}\right) - H_y^{n+\frac{1}{2}}\!\left(i,k-\frac{1}{2}\right) }{\Delta z}. \]
Reúne los términos \(E_x^{n+1}\) en el lado izquierdo:
\[ \left( \frac{\sigma}{2} + \frac{\epsilon}{\Delta t} \right) E_x^{n+1}(i,k) = - \frac{ H_y^{n+\frac{1}{2}}\!\left(i,k+\frac{1}{2}\right) - H_y^{n+\frac{1}{2}}\!\left(i,k-\frac{1}{2}\right) }{\Delta z} + \left( \frac{\epsilon}{\Delta t} - \frac{\sigma}{2} \right) E_x^n(i,k). \]
Por lo tanto.
\[ E_x^{n+1}(i,k) = \frac{ \left( \dfrac{\epsilon}{\Delta t} - \dfrac{\sigma}{2} \right) E_x^n(i,k) - \dfrac{ H_y^{n+\frac{1}{2}}\!\left(i,k+\frac{1}{2}\right) - H_y^{n+\frac{1}{2}}\!\left(i,k-\frac{1}{2}\right) }{\Delta z} }{ \dfrac{\epsilon}{\Delta t} + \dfrac{\sigma}{2} }. \]
Multiplicar el numerador y el denominador por \(2\Delta t\) da como resultado la forma de coeficiente que ya conocemos:
\[ E_x^{n+1}(i,k) = \frac{ 2\epsilon - \sigma \Delta t }{ 2\epsilon + \sigma \Delta t } E_x^n(i,k) - \frac{ 2 \Delta t }{ 2\epsilon + \sigma \Delta t } \frac{ H_y^{n+\frac{1}{2}}\!\left(i,k+\frac{1}{2}\right) - H_y^{n+\frac{1}{2}}\!\left(i,k-\frac{1}{2}\right) }{\Delta z}. \]
6.2. Actualizar la ecuación para \(E_y\).
Comienza desde
\[ \sigma E_y + \epsilon \frac{\partial E_y}{\partial t} = - \frac{\partial H_z}{\partial x} + \frac{\partial H_x}{\partial z}. \]
Discretizar ambos términos de la derivada en la malla de Yee:
\[ \sigma \frac{E_y^{n+1}(i,k) + E_y^n(i,k)}{2} + \epsilon \frac{E_y^{n+1}(i,k) - E_y^n(i,k)}{\Delta t} = - \frac{ H_z^{n+\frac{1}{2}}\!\left(i+\frac{1}{2},k\right) - H_z^{n+\frac{1}{2}}\!\left(i-\frac{1}{2},k\right) }{\Delta x} + \frac{ H_x^{n+\frac{1}{2}}\!\left(i,k+\frac{1}{2}\right) - H_x^{n+\frac{1}{2}}\!\left(i,k-\frac{1}{2}\right) }{\Delta z}. \]
Reorganizar permite obtener...
\[ E_y^{n+1}(i,k) = \frac{ 2\epsilon - \sigma \Delta t }{ 2\epsilon + \sigma \Delta t } E_y^n(i,k) + \frac{ 2 \Delta t }{ 2\epsilon + \sigma \Delta t } \left[ - \frac{ H_z^{n+\frac{1}{2}}\!\left(i+\frac{1}{2},k\right) - H_z^{n+\frac{1}{2}}\!\left(i-\frac{1}{2},k\right) }{\Delta x} + \frac{ H_x^{n+\frac{1}{2}}\!\left(i,k+\frac{1}{2}\right) - H_x^{n+\frac{1}{2}}\!\left(i,k-\frac{1}{2}\right) }{\Delta z} \right]. \]
Este es el único componente del campo eléctrico en la reducción bidimensional que contiene contribuciones de dos derivadas distintas del campo magnético, y por lo tanto es el más fácil de calcular incorrectamente en el código. En particular, el signo relativo entre los términos de las derivadas \(x\) y \(z\) se deriva directamente de la expansión del rotacional y no debe modificarse.
6.3. Actualizar la ecuación para \(E_z\).
Comienza desde
\[ \sigma E_z + \epsilon \frac{\partial E_z}{\partial t} = \frac{\partial H_y}{\partial x}. \]
La discretización permite obtener
\[ \sigma \frac{E_z^{n+1}(i,k) + E_z^n(i,k)}{2} + \epsilon \frac{E_z^{n+1}(i,k) - E_z^n(i,k)}{\Delta t} = \frac{ H_y^{n+\frac{1}{2}}\!\left(i+\frac{1}{2},k\right) - H_y^{n+\frac{1}{2}}\!\left(i-\frac{1}{2},k\right) }{\Delta x}. \]
Por lo tanto.
\[ E_z^{n+1}(i,k) = \frac{ 2\epsilon - \sigma \Delta t }{ 2\epsilon + \sigma \Delta t } E_z^n(i,k) + \frac{ 2 \Delta t }{ 2\epsilon + \sigma \Delta t } \frac{ H_y^{n+\frac{1}{2}}\!\left(i+\frac{1}{2},k\right) - H_y^{n+\frac{1}{2}}\!\left(i-\frac{1}{2},k\right) }{\Delta x}. \]
En este punto, ya se han obtenido las ecuaciones discretas para la actualización del campo eléctrico, presentadas en su forma estándar con amortiguación. El siguiente paso consiste en obtener las ecuaciones correspondientes para el campo magnético, utilizando niveles de tiempo semi-enteros.
7. Obtención de las ecuaciones para la actualización del campo magnético.
El campo magnético se actualiza en intervalos de tiempo que son múltiplos impares de la mitad. Con \(\sigma_m = 0\), las ecuaciones del campo magnético no contienen ningún término de amortiguación, por lo que las fórmulas de actualización adoptan una forma más sencilla y directa. Esta parte de la derivación es algebraicamente más corta que el caso del campo eléctrico, ya que no hay ningún término de conductividad trapezoidal que deba incluirse.
7.1 Actualizar la ecuación para \(H_x\).
Comienza desde
\[ \mu \frac{\partial H_x}{\partial t} = \frac{\partial E_y}{\partial z}. \]
Al discretizar en el tiempo entre los puntos medios \(n-\tfrac{1}{2}\) y \(n+\tfrac{1}{2}\) y en el espacio utilizando una diferencia centrada en la dirección \(z\), se obtiene:
\[ \mu \frac{ H_x^{n+\frac{1}{2}}(i,k+\frac{1}{2}) - H_x^{n-\frac{1}{2}}(i,k+\frac{1}{2}) }{\Delta t} = \frac{ E_y^n(i,k+1) - E_y^n(i,k) }{\Delta z}. \]
En la notación de índice medio reducido, que se utiliza con frecuencia en las demostraciones, el mismo resultado se puede expresar de la siguiente manera:
\[ \mu \frac{ H_x^{n+\frac{1}{2}} - H_x^{n-\frac{1}{2}} }{\Delta t} = \frac{ E_y^n\!\left(k+\frac{1}{2}\right) - E_y^n\!\left(k-\frac{1}{2}\right) }{\Delta z}. \]
Por lo tanto.
\[ H_x^{n+\frac{1}{2}} = H_x^{n-\frac{1}{2}} + \frac{\Delta t}{\mu} \frac{ E_y^n\!\left(k+\frac{1}{2}\right) - E_y^n\!\left(k-\frac{1}{2}\right) }{\Delta z}. \]
7.2. Actualizar la ecuación para \(H_y\).
Comienza desde
\[ \mu \frac{\partial H_y}{\partial t} = \frac{\partial E_z}{\partial x} - \frac{\partial E_x}{\partial z}. \]
La discretización permite obtener
\[ \mu \frac{ H_y^{n+\frac{1}{2}}(i+\frac{1}{2},k+\frac{1}{2}) - H_y^{n-\frac{1}{2}}(i+\frac{1}{2},k+\frac{1}{2}) }{\Delta t} = \frac{ E_z^n(i+1,k+\frac{1}{2}) - E_z^n(i,k+\frac{1}{2}) }{\Delta x} - \frac{ E_x^n(i+\frac{1}{2},k+1) - E_x^n(i+\frac{1}{2},k) }{\Delta z}. \]
De nuevo, utilizando la notación simplificada del índice parcial:
\[ \mu \frac{ H_y^{n+\frac{1}{2}} - H_y^{n-\frac{1}{2}} }{\Delta t} = \frac{ E_z^n\!\left(i+\frac{1}{2}\right) - E_z^n\!\left(i-\frac{1}{2}\right) }{\Delta x} - \frac{ E_x^n\!\left(k+\frac{1}{2}\right) - E_x^n\!\left(k-\frac{1}{2}\right) }{\Delta z}. \]
Por lo tanto.
\[ H_y^{n+\frac{1}{2}} = H_y^{n-\frac{1}{2}} + \frac{\Delta t}{\mu} \left[ \frac{ E_z^n\!\left(i+\frac{1}{2}\right) - E_z^n\!\left(i-\frac{1}{2}\right) }{\Delta x} - \frac{ E_x^n\!\left(k+\frac{1}{2}\right) - E_x^n\!\left(k-\frac{1}{2}\right) }{\Delta z} \right]. \]
7.3. Actualización de la ecuación para \(H_z\).
Comienza desde
\[ \mu \frac{\partial H_z}{\partial t} = - \frac{\partial E_y}{\partial x}. \]
La discretización permite obtener
\[ \mu \frac{ H_z^{n+\frac{1}{2}}(i+\frac{1}{2},k) - H_z^{n-\frac{1}{2}}(i+\frac{1}{2},k) }{\Delta t} = - \frac{ E_y^n(i+1,k) - E_y^n(i,k) }{\Delta x}. \]
De manera equivalente, en una notación más concisa:
\[ \mu \frac{ H_z^{n+\frac{1}{2}} - H_z^{n-\frac{1}{2}} }{\Delta t} = - \frac{ E_y^n\!\left(i+\frac{1}{2}\right) - E_y^n\!\left(i-\frac{1}{2}\right) }{\Delta x}. \]
Por lo tanto.
\[ H_z^{n+\frac{1}{2}} = H_z^{n-\frac{1}{2}} - \frac{\Delta t}{\mu} \frac{ E_y^n\!\left(i+\frac{1}{2}\right) - E_y^n\!\left(i-\frac{1}{2}\right) }{\Delta x}. \]
Una vez realizadas estas sustituciones, se pueden aplicar directamente las ecuaciones estándar de actualización del método FDTD.
8. Ecuaciones finales de actualización.
Una vez recopilados los resultados, las ecuaciones de actualización bidimensionales de FDTD son:
\[ E_x^{n+1}(i,k) = \frac{ 2\epsilon - \sigma \Delta t }{ 2\epsilon + \sigma \Delta t } E_x^n(i,k) - \frac{ 2 \Delta t }{ 2\epsilon + \sigma \Delta t } \frac{ H_y^{n+\frac{1}{2}}\!\left(i,k+\frac{1}{2}\right) - H_y^{n+\frac{1}{2}}\!\left(i,k-\frac{1}{2}\right) }{\Delta z}, \]
\[ E_y^{n+1}(i,k) = \frac{ 2\epsilon - \sigma \Delta t }{ 2\epsilon + \sigma \Delta t } E_y^n(i,k) + \frac{ 2 \Delta t }{ 2\epsilon + \sigma \Delta t } \left[ - \frac{ H_z^{n+\frac{1}{2}}\!\left(i+\frac{1}{2},k\right) - H_z^{n+\frac{1}{2}}\!\left(i-\frac{1}{2},k\right) }{\Delta x} + \frac{ H_x^{n+\frac{1}{2}}\!\left(i,k+\frac{1}{2}\right) - H_x^{n+\frac{1}{2}}\!\left(i,k-\frac{1}{2}\right) }{\Delta z} \right], \]
\[ E_z^{n+1}(i,k) = \frac{ 2\epsilon - \sigma \Delta t }{ 2\epsilon + \sigma \Delta t } E_z^n(i,k) + \frac{ 2 \Delta t }{ 2\epsilon + \sigma \Delta t } \frac{ H_y^{n+\frac{1}{2}}\!\left(i+\frac{1}{2},k\right) - H_y^{n+\frac{1}{2}}\!\left(i-\frac{1}{2},k\right) }{\Delta x}, \]
\[ H_x^{n+\frac{1}{2}}(i,k+\tfrac{1}{2}) = H_x^{n-\frac{1}{2}}(i,k+\tfrac{1}{2}) + \frac{\Delta t}{\mu} \frac{ E_y^n(i,k+1) - E_y^n(i,k) }{\Delta z}, \]
\[ H_y^{n+\frac{1}{2}}(i+\tfrac{1}{2},k+\tfrac{1}{2}) = H_y^{n-\frac{1}{2}}(i+\tfrac{1}{2},k+\tfrac{1}{2}) + \frac{\Delta t}{\mu} \left[ \frac{ E_z^n(i+1,k+\tfrac{1}{2}) - E_z^n(i,k+\tfrac{1}{2}) }{\Delta x} - \frac{ E_x^n(i+\tfrac{1}{2},k+1) - E_x^n(i+\tfrac{1}{2},k) }{\Delta z} \right], \]
\[ H_z^{n+\frac{1}{2}}(i+\tfrac{1}{2},k) = H_z^{n-\frac{1}{2}}(i+\tfrac{1}{2},k) - \frac{\Delta t}{\mu} \frac{ E_y^n(i+1,k) - E_y^n(i,k) }{\Delta x}. \]
Tal como se describe aquí, el significado del esquema iterativo es evidente. El campo eléctrico en el instante \(n+1\) se obtiene a partir del campo magnético en el instante \(n+\tfrac{1}{2}\), mientras que el campo magnético en el instante \(n+\tfrac{1}{2}\) se obtiene a partir del campo eléctrico en el instante \(n\). Por lo tanto, los campos evolucionan en una secuencia alternada en el tiempo.
9. Observaciones sobre la implementación.
En un código práctico, los coeficientes que multiplican los valores antiguos del campo y los términos de rotacional suelen calcularse y almacenarse de antemano, especialmente cuando \(\epsilon\), \(\mu\) o \(\sigma\) varían de una celda a otra. Por ejemplo, la actualización del campo eléctrico suele adoptar la siguiente forma:
\[ E^{n+1} = C_a E^n + C_b (\nabla \times H), \]
donde los coeficientes locales \(C_a\) y \(C_b\) codifican los parámetros del material y el intervalo de tiempo. Esto simplifica la implementación y evita tener que recalcular los mismos factores en cada intervalo de tiempo.
También es importante destacar que la derivación anterior solo proporciona el esquema básico de actualización de FDTD. Para realizar una simulación útil, también se necesitan fuentes, condiciones de contorno, capas absorbentes, como las capas de adaptación perfecta (PML), y un paso de tiempo que cumpla con la condición de estabilidad de Courant. Ninguna de estas adiciones altera la estructura fundamental del método de «salto de rana» que se deriva aquí, pero cada una de ellas modifica la forma en que se aplican las ecuaciones de actualización en ubicaciones específicas de la malla.
Por último, aunque la notación anterior se ha redactado con suficiente cuidado para que sea relevante para su implementación, la ubicación espacial exacta de cada componente depende de la convención de cuadrícula utilizada en el código. El requisito fundamental no es la tipografía en sí, sino la coherencia: el operador de rizo, el escalonamiento y los signos deben corresponder todos a la misma disposición de celdas de Yee en todo el código.