Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
4
1. Propagación en el espacio libre
En este capítulo hacemos una presentación sencilla del MDFDT discutiendo la
formulación para simular la propagación del campo electromagnético en un medio
homogéneo. Ciertamente, la propagación en un medio homogéneo o en el vacío no es un
problema que implique un gran interés físico. Sin embargo, el objetivo es la ilustración
del algoritmo.
1.1 Las ecuaciones de Maxwell en forma discreta
EL MDFDT analiza el problema de la propagación electromagnética en pequeñas
particiones espaciales. En estas celdas los campos eléctricos y magnéticos están
alternadamente distribuidos. Las celdas o nodos también estas intercaladas en el tiempo.
Para resolver este conjunto de ecuaciones, se desarrolla un proceso iterativo en el tiempo.
Comenzamos con las ecuaciones de Maxwell en el sistema (cgs).5
JxDxHct
tc
t π4),(1),( +∂
∂=×∇ (1.1-1)
t
tc
t∂
∂−=×∇
),(1),( xBxE (1.1-2)
πρ4),( =⋅∇ txD (1.1-3)
0),( =⋅∇ txB (1.1-4)
Donde B corresponde a la inducción magnética, E corresponde al campo
eléctrico, H corresponde al campo magnético, D es el desplazamiento eléctrico, ρ es la
densidad de carga, J es la densidad de corriente eléctrica y c corresponde a la velocidad
de la luz en el vacío. Debemos hacer notar que para medios isotrópicos y lineales las
relaciones constitutivas se escriben:
ED ε=
HB μ=
EJ .σ=
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
5
Siendo ε la constante dieléctrica del material, μ la permeabilidad magnética del material y
σ la conductividad.
Para medios dispersivos
)()()( ωωεω ED =
)()()( ωωμω HB =
)()()( ωωσω EJ =
Para ondas que se propagan en medios homogéneos y en dieléctricos y en donde
no existen fuentes de carga ni de corrientes (ρ = J = 0), las ecuaciones (1.1-1) y (1.1-2)
toman la forma:
tt
ct
∂∂
=×∇),(),( xExH ε
(1.2-1)
t
tc
t∂
∂−=×∇
),(),( xHxE μ (1.2-2)
Si escogemos la dirección x para la polarización del campo eléctrico y la
dirección z para la dirección de la propagación, entonces el campo magnético tiene la
polarización y. Las ecuaciones rotacionales vectoriales (1.2-1) y (1.2-2) son ahora las
ecuaciones escalares,
),(),( tzHz
ctzEt yx ∂
∂−=
∂∂
ε (1.3-1)
),(),( tzEz
ctzHt xy ∂
∂−=
∂∂
μ (1.3-2)
La primera de estas ecuaciones es también llamada de Maxwell-Ampere y la
segunda se conoce como la ecuación de Faraday. Estas ecuaciones son ecuaciones
puntuales, es decir, que son válidas para cada valor continuo de las coordenadas espacial
y temporal. Para lograr esto tenemos que trabajar en un espacio discreto en lugar de
hacerlo en un espacio continuo.
Tomemos la definición de derivada en la variante de diferencias centrales un
punto zo,6
z
zzfzzfdz
zdf oo
oz ΔΔ−−Δ+
=→Δ
)2/()2/(lim
)( 0 (1.3-3)
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
6
f(z - z )
z0
f(z )0
z 0
0 2
z
f(z + z )0 2
f (z )
(z - z )20 0(z + z )
2
f '(z )0
Esta ecuación es valida para todos los puntos del espacio z. En esta ecuación, para un
valor finito de Δz, la ecuación diferencial se transforma en una ecuación de diferencias
finitas. La derivada esta definida por el valor de la función en puntos discretos contiguos
a zo.
Figura 1.1. Representación gráfica de derivada en diferencias centrales.
Es necesario discretizar las ecuaciones de Maxwell. La estrategia consiste en
definir un punto “z” por medio de un valor discreto “z = kΔz” y de igual forma la
coordenada temporal “t” toma valores discretos de la forma “t = n Δt”. Para ilustrar estas
ideas obsérvese la figura en donde los valores accesibles de z y t están determinados por
los saltos de rana discretos que se muestran en la Fig. 1.2, cada salto de rana simula un
paso espacial o temporal según sea el caso.
Figura 1.2. Representación de la utilización de valores discretos de z y t.
kΔz (k+1)Δz z
nΔt (n+1)Δt t
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
7
Con estas consideraciones, la forma discreta de la derivada espacial para el campo
eléctrico (E) es:
ztnzkEtnzkE
tnzkEz
xxx Δ
ΔΔ−−ΔΔ+≅ΔΔ
∂∂ ],)2/1[(],)2/1[(
),(
Similarmente, los valores discretos de la coordenada temporal están indexados por el
entero ‘n’ por medio de la relación t = nΔt y la forma discreta temporal queda
ttnzkEtnzkE
tnzkEt
xxx Δ
Δ−Δ−Δ+Δ≅ΔΔ
∂∂ ])2/1(,[])2/1(,[
),(
Vamos a usar la convención de escribir en el superíndice el valor discreto de la
variable temporal y entre paréntesis el valor discreto de la variable espacial, es decir,
escribimos
)(),(),( kEtnzkEtzE nxxx =ΔΔ≅
Relaciones similares pueden obtenerse para el campo Hy. Al implementar las
ecuaciones en diferencias finitas para Ex(z,t) y Hy(z,t) en la ec. (1.1-7) se obtiene la forma
discreta de la ecuación de Maxwell-Ampere:
z
kHkHct
kEkE ny
ny
nx
nx
Δ
−−+−=
Δ− −+ )2/1()2/1()()( 2/12/1
ε (1.4-1)
Al igual que hemos hecho la consideración de los campos en un punto (z,t), ahora
vamos a considerar la siguiente ecuación rotacional (1.1-8) en el punto
)2/,2/( ttzz Δ+Δ+ , así la ecuación de Faraday nos queda:
z
kEkEct
kHkH nx
nx
ny
ny
Δ−+
−=Δ
+−+ +++ )()1()2/1()2/1( 2/12/11
μ (1.4-2)
Al considerar a los campos en coordenadas diferentes, es decir (z,t) y
)2/,2/( ttzz Δ+Δ+ , logramos describir una situación de evolución espacio temporal.
Estas ecuaciones pueden expresarse de la forma:
[ ])2/1()2/1()()( 2/12/1 −−+ΔΔ
−= −+ kHkHztckEkE n
yny
nx
nx ε
(1.5-1)
[ ])()1()2/1()2/1( 2/12/11 kEkEztckHkH n
xnx
ny
ny
+++ −+ΔΔ
−+=+μ
(1.5-2)
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
8
Estas ecuaciones están intercaladas en ambas coordenadas: espacio y tiempo. Esta forma
de escribir las ecuaciones de Maxwell tiene la ventaja que asemeja la forma usual en que
visualizamos la propagación de un campo electromagnético: un campo eléctrico induce
a un campo magnético, el cual induce un campo eléctrico, y así al infinito. Este
mecanismo es esquematizado en la Fig. 1.3. En esta figura, las flechas indican la forma
en como el nuevo valor del campo eléctrico proviene del campo eléctrico en ese punto un
paso temporal anterior más la contribución de los campos Hy en puntos contiguos en el
espacio un medio tiempo atrás. Este es el paradigma fundamental del MDFDT.
Hn(k+3/2)Hn(k-1/2)Hn(k-3/2) Hn(k+1/2)
En-1/2(k-2) En-1/2(k-1) En-1/2(k+2)En-1/2(k+1)
(n+2)Δt
(n+3/2)Δt
(n+1)Δt
(n+1/2)Δt
nΔt
kΔz-2Δz kΔz-Δz kΔz+2Δz kΔz+Δz
tiem
po (t
)
distancia (z)
E H
kΔz
(n-1/2)ΔtEn-1/2(k)
En+1/2(k-2) En+1/2(k-1) En+1/2(k+2)En+1/2(k+1)En+1/2(k)
Hn+1(k+3/2)Hn+1(k-1/2)Hn+1(k-3/2) Hn+1(k+1/2)
En+3/2(k-2) En+3/2(k-1) En+3/2(k+2)En+3/2(k+1)En+3/2(k)
Figura 1.3. Representación de la evolución temporal de un campo electromagnético.
En este punto, es necesario determinar las cantidades que afectan a los corchetes
en las ecs. (1.5-1) y (1.5-2). En el vacío, la función dieléctrica y la permeabilidad
magnética son iguales a 1; ε = μ =1. Una onda electromagnética necesita un mínimo de
tiempo Δt =Δz/c para propagarse dentro de la malla discreta espacial Δz. Si tomamos un
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
9
paso temporal mas grande, la onda ya habría pasado la distancia Δz y ya no estaría en
nuestra malla de simulación. Por esta razón es necesario escoger un paso temporal de la
forma: 1,7
czt Δ
≤Δ (1.6-1)
En este trabajo vamos a utilizar la convención de escoger un paso temporal
czt
2Δ
=Δ (1.7-1)
Para la propagación en el vacío, el paso espacial lo escogemos tomando 10
puntos de la longitud de onda de la luz en el vacío λ0,
100λ
=Δz . (1.8-1)
1.2 Implementación Computacional
El primer campo electromagnético que vamos a simular es un pulso gaussiano que
tiene la forma:
⎪⎪
⎩
⎪⎪
⎨
⎧
>
<=
⎟⎟⎠
⎞⎜⎜⎝
⎛ −−
0
00)(
20
21
te
ttE
t
ttσ
Esta es una función temporal que asemeja un pulso de luz y esta ilustrada en la Fig.
1.4. Tiene la característica de que su amplitud máxima es E0. Solo toma valores
diferentes de cero en la vecindad de t0, tiempo al cual esta centrada y tiene una duración
media σ t.
La implementación de las ecs. (1.5-1) y (1.5-2) en un algoritmo computacional
requiere dos vectores unidimensionales, uno para Ex y otro para Hy, para describir la
coordenada espacial. En cada tiempo ‘n’ los campos son calculados en Matlab mediante
las líneas de programa:
ex(k) = ex(k) + 0.5*( hy(k-1) – hy(k) ) (1.9-1)
hy(k) = hy(k) + 0.5*( ex(k) – ex(k+1) ) (1.9-2)
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
10
0 1 2 3 4 5 6 7 8 9 100,0
0,2
0,4
0,6
0,8
1,0
t0=5τ
σt
C
ampo
ele
ctric
o (1
/E0)
tiempo (t/τ)
Figura 1.4 Ilustración de un pulso gaussiano.
En estas relaciones, ya no existe el superíndice. El tiempo es una variable
implícita en el MDFDT. En la primera de estas ecuaciones, el nuevo valor de ex(k) a la
izquierda de la relación [para el tiempo (n+1/2)] proviene del antiguo valor ex(k) [al
tiempo(n-1/2)] que la máquina ya ha calculado y tenia en la memoria más los campos
magnéticos hy(k-1)y hy(k) [propios del tiempo n]. Para el vector hy los índices
se han redondeado de los índices k+1/2 y k-1/2 a k y k-1 para poder especificar la
posición de los campos en el arreglo computacional.
En el apéndice 1 mostramos el programa ‘simple.m’ que es la primera
implementación del algoritmo MDFDT en Matlab. Ya que el MDFDT es un proceso
iterativo, puede ser presentado en forma concisa en un diagrama de flujo.
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
11
Como se aprecia en el diagrama el programa comienza inicializando constantes y
parámetros. Los campos eléctricos y magnéticos son inicializados a cero. El cálculo
comienza obteniendo los campos eléctricos en cada nodo usando la relación (1.6-1). A
continuación, se introduce el valor de la fuente de campo eléctrico. El siguiente paso es
calcular los campos magnéticos. Finalmente, el programa grafica e itera el procedimiento
al siguiente valor temporal.
En la Fig. 1.5 mostramos el resultado del programa para 100 iteraciones
temporales. El comportamiento del pulso en el espacio y en el tiempo puede ser
visualizado en la Fig. 1.6. A medida que el pulso se propaga en el tiempo, se extiende
hacia fuera del centro de la malla espacial.
Definición de parámetros e inicialización
de los campos
Cálculo del campo eléctrico
Inserción de la fuente
Cálculo del campo magnético
Graficación
Fín
Inicio
¿Fin de pasos
temporales?
No
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
12
El valor de la longitud de onda utilizado en las simulaciones para este capítulo es
del orden de 632.8x10-9 m.
El programa detiene la simulación antes de tocar el límite de la malla de
simulación. En este método, para calcular el campo E en un punto, es necesario conocer
los valores H que rodean dicho punto. Sin embargo, en el borde de la malla de simulación
no conocemos el punto correspondiente que esta fuera de la malla. De esta forma, no es
posible calcular el campo E en los valores límites de la malla.
Si el programa continuase, es decir, si hacemos el ejercicio de modificar el
número de pasos temporales a N_t=200, veremos en la evolución temporal los valores
falsos que toma el campo en el borde de la malla. El efecto neto es una reflexión hacia el
interior de la malla. Esto es algo que es necesario evitar ya que esto no es físicamente
aceptable. Por ello introducimos las condiciones de frontera absorbentes.
Figura 1.5. Propagación de un pulso gaussiano en el espacio libre.
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
13
Figura 1.6. Variación de un pulso gaussiano en el espacio y en el tiempo
1.3 Condiciones de frontera absorbentes
Las condiciones de frontera absorbentes son necesarias para simular
adecuadamente la propagación de los campos electromagnéticos una vez que estos han
llegado a la frontera de nuestra malla de simulación. En el algoritmo del MDFDT, los
valores de los campos se determinan mediante un promedio de los campos en los puntos
vecinos. El problema consiste en que en la frontera de la malla de simulación, este
promedio no se puede dar porque no conocemos el valor del campo fuera de la malla. De
esta forma, si la simulación continua, los campos toman valores no válidos.
Físicamente, esperamos que los campos se propaguen hacia fuera de la malla de
simulación, ya que consideramos que fuera de nuestra malla no existen fuentes. La
distancia que la onda viaja en un paso temporal esta determinado por la relación (1.7-1).
distancia = c0Δt=c0(Δz/2c0)=Δz/2 (1.7-1)
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
14
Figura 1.7 Pulso gaussiano después de 140 pasos temporales incluyendo
condiciones absorbentes de frontera.
Esta relación implica que al campo le toma un par de pasos temporales viajar un
paso espacial. La forma de expresar esta condición es
)2()1( 2−= nx
nx EE
Con esta relación estamos asignando el valor del campo en la frontera en lugar de
calcularlo. Esta condición es fácil de implementar. Solamente es necesario guardar el
valor de Ex(2) un par de pasos temporales y luego asignarlos a Ex(1). Las condiciones de
frontera han sido implementadas en el programa ‘simple_ABC.m’ que se encuentra en el
apéndice 2. En la figura 1.7 mostramos la propagación de un pulso gaussiano. En esta
ocasión, el pulso pasa la frontera sin reflexión, que es lo esperado que ocurra.
1.4 Aplicaciones
Una de las ventajas más significativas del MDFDT es que es flexible para realizar
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
15
simulaciones tanto de ondas como de pulsos que varían en el tiempo. También es posible
simular a varias frecuencias, por ejemplo una onda sinusoidal o un pulso que contenga
varias frecuencias. A continuación veremos algunos ejemplos que ilustran la flexibilidad
del algoritmo.
1.4.1 Superposición de dos pulsos gaussianos
El primer ejemplo de modificación del programa es obtenido al analizar la
superposición de dos pulsos gaussianos. Uno de los pulsos esta ubicado en k=30 y el
segundo esta ubicado en k=70. En la Fig. 1.8 ilustramos la superposición de ambos
pulsos después de 90 pasos temporales, en el momento que se detiene la simulación los
máximos de cada pulso coinciden, aquí la amplitud total es la suma de las amplitudes
individuales, posterior a su encuentro cada pulso continua su trayecto.
Para lograr esta simulación, es necesario escribir una nueva función llamada
‘pulso_2gauss.m’ la cual mostramos en el apéndice 3. Asimismo, hay que modificar una
línea del programa ‘simple_ABC.m’. La línea 28 donde viene la instrucción
ex=pulso(ex,n,t,tao,t0,spread);
Debe ser cambiada por
ex=pulso_2gauss(ex,n,t,tao,t0,spread);
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
16
Figura 1.8 Superposición de dos pulsos gaussianos
Con esta variación estamos llamando a la nueva función ‘pulso_2gauss.m’ en vez
de a la función ‘pulso.m’. Cabe destacar que tanto el programa, como la función que
estemos llamando siempre deben de estar en el mismo subdirectorio.
1.4.2 Superposición de dos ondas sinusoidales
Este ejemplo requiere mas cambios sobre el programa ‘simple_ABC’, Ahora
vamos a requerir de una malla temporal mas grande, N_z=400 y las iteraciones
temporales van a ser N_t=350. La onda sinusoidal que vamos a simular es de la forma:
⎟⎠⎞
⎜⎝⎛= ttE
τπ21.0sin)(
Esta fuente va a ser ubicada en k=20 y k=380. Este cambio esta contenido en la
función ‘pulso_2sin.m’, la cual se llama dentro del programa modificando la línea 28
donde viene la instrucción
ex=pulso(ex,n,t,tao,t0,spread);
debe ser cambiada por
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
17
ex=pulso_2sin(ex,n,t,tao,t0,spread);
En la Fig. 1.9 mostramos la superposición de dos ondas después de 405 pasos
temporales, lo único que difiere de la anterior simulación es el tipo de fuente que se está
utilizando, la simulación se detiene cuando los máximos de las amplitudes en las ondas
coinciden.
Figura 1.9 Superposición de dos ondas sinusoidales.
1.4.3 Evolución de un grupo de ondas
Cuando dos ondas se encuentran presentes en el mismo medio estas pueden
combinarse constructivamente o destructivamente dependiendo de la coincidencia de los
máximos.
Para ilustrar este caso supongamos dos ondas sinusoidales de la forma:
⎟⎟⎠
⎞⎜⎜⎝
⎛= ttE
τπ2
3.0sin)(1
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
18
⎟⎟⎠
⎞⎜⎜⎝
⎛= ttE
τπ2
28.0sin)(2
El grupo de ondas formado por los campos anteriores obedece a la función:
⎟⎠⎞
⎜⎝⎛+⎟
⎠⎞
⎜⎝⎛= tttE
τπ
τπ 228.0sin23.0sin)(
Similarmente, la línea 28 donde viene la instrucción
ex=pulso(ex,n,t,tao,t0,spread);
debe ser cambiada por
ex=pulso_2group (ex,n,t,tao,t0,spread);
En este ejemplo consideramos una malla de tamaño N_z = 500 y en la fig. 2.10
mostramos la situación después de 1000 pasos temporales.
Figura 1.10 Evolución de la suma de dos ondas sinuosidades ligeramente desfasadas.
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
19
0 1 2 3
x 10-4
-1
-0.5
0
0.5
1
z (metros)
E(z
,t)
0 1 2 3
x 10-4
-1
-0.5
0
0.5
1
z (metros)
H(z
,t)
1.4.4 Dirección de propagación.
La dirección de propagación de la luz esta dictada por la regla de la mano derecha
según lo indica el vector de Poynting:
HES ×=π4c .
Dependiendo del signo que acompañe a los campos esto determinará la dirección del
pulso. Sea entonces la fuente del campo eléctrico para t=0 de la forma 22/
21
sin)(⎟⎟⎠
⎞⎜⎜⎝
⎛
Δ
−−
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛= z
zz
fx
f
ezz
zE π
22/21
sin)(⎟⎟⎠
⎞⎜⎜⎝
⎛
Δ
−−
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛= z
zz
fy
f
ezz
zH π
Cuando ambos campos son tomados positivos, estos se propagan en la dirección positiva,
como lo muestra la fig. 1.11.
Figura 1.11 Ilustración de una onda electromagnética propagándose a la derecha.
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
20
Cuando alguno de los campos es negativo, tomemos como ejemplo al campo eléctrico
negativo entonces la onda correspondiente se propaga en la dirección negativa como lo
ilustra la fig. 1.12.
Figura 1.12 Ilustración de una onda electromagnética propagándose a la izquierda.
1.4.5 Simulación de una cavidad
Como un último ejemplo, vamos a simular la oscilación dentro de una cavidad.
En este caso, el campo es cero en la frontera de la malla de simulación, lo cual
físicamente representa una cavidad de aire en medio de paredes perfectamente
conductoras. La función que describe el campo es
⎟⎟⎠
⎞⎜⎜⎝
⎛= z
zzE
fx
πsin)(
donde zf define el mayor valor asignado para la coordenada espacial.
0 1 2 3
x 10-4
-1
-0.5
0
0.5
1
z (metros)
E(z
,t)
0 1 2 3
x 10-4
-1
-0.5
0
0.5
1
z (metros)
H(z
,t)
Aplicación del Método de Diferencias Finitas en el Dominio del Tiempo
a la simulación del Campo Electromagnético.
21
Este campo esta definido solamente en el tiempo inicial t=0, a continuación, el algoritmo
se encarga de realizar la simulación de la onda que oscila en el tiempo. Esta cavidad se
ilustra en la fig. 1.13.
Figura 1.13 Ilustración de una oscilación dentro de una cavidad.
Top Related