View
3
Download
2
Category
Preview:
DESCRIPTION
MATLAB
Citation preview
Relatório: Fotónica 1
1. Resolução numérica de equações diferenciais ordinárias Na simulação numérica da modulação directa da corrente de injecção em lasers
semicondutores há que resolver as equações das taxas. Considera-se, aqui, o caso dos lasers
semicondutores monomodais. Estas equações regulam a evolução temporal, na zona activa do
laser semicondutor, das populações totais de electrões N (envolvidos nos processos de
geração e recombinação) e de fotões S , tendo-se
equações das taxasst sp a
p t
d S R R Rdtd N R Rdt
′= Γ + −→
= −
em que 1Γ < representa o factor de confinamento óptico (o modo de oscilação laser não fica
exactamente confinado à zona activa). As várias taxas intervenientes têm unidades ( )1s− e
são: (i) a taxa de bombeamento pR ; (ii) a taxa total de recombinação (de pares electrão-
lacuna) tR ; (iii) a taxa efectiva líquida de emissão estimulada stR ; (iv) a taxa contributiva de
emissão espontânea spR′ ; (v) a taxa efectiva de aniquilação de fotões aR (relacionada com o
tempo de vida médio pτ dos fotões na cavidade laser).
Shun Lien Chuang, Physics of Optoelectronic Devices (New York: Wiley, 1995),
Chapter 11 (pp. 487-507).
Tem-se
2 Carlos R. Paiva
, , , ,p t st st sp sp sp ap
I NR R R R R G S R R Rq
βτ
′= = + = = =
em que I representa a corrente de injecção (sendo q o valor da carga do electrão,
considerada aqui como positiva), G a taxa elementar líquida de emissão estimulada e spβ o
coeficiente de emissão espontânea contributiva (i.e., a pequena parcela da taxa de emissão
espontânea total spR que contribui para a emissão estimulada). A taxa R , que integra a taxa
total de recombinação tR , representa a taxa de aniquilação de electrões e é dada por
c
NRτ
= .
Em geral cτ , que é o tempo de vida médio dos electrões, depende de N . Para simplificar a
análise considera-se, frequentemente, que cτ é uma constante. A taxa spR′ relaciona-se com
G através da expressão
sp sp sp spR R n Gβ′ = =
em que spn é o chamado coeficiente de inversão da população. As equações das taxas
assumem, portanto, a forma explícita
equações das taxassp
p
c
d S SG S n Gdt
d N I NG Sdt q
τ
τ
= Γ + −
→= − −
e, para serem resolvidas, necessitam que se adopte um modelo que relacione G com N e S .
Um modelo possível (não linear) é
( )1
N tG N NG
Sε−
=+
Relatório: Fotónica 3
onde NG , tN e ε são constantes que caracterizam este modelo específico. A resolução
numérica das equações das taxas enquadra-se, portanto, no quadro mais geral da resolução
numérica de equações diferenciais ordinárias.
Apresenta-se, a seguir, um caso concreto de um sistema de duas equações diferenciais
ordinárias que pertence à mesma classe das equações das taxas. Considere-se, com efeito, que
( ) ( )
( ) ( )
( ) ( )( ) ( )
1 111 2
1 2
2 222 1
1 2
11
11
B p zd p A p zd z p z p z
B p zd p A p zd z p z p z
⎧= −⎡ ⎤⎪ ⎣ ⎦ + +⎪
⎨⎪ = −⎡ ⎤⎣ ⎦⎪ + +⎩
onde 1 1 2 2, , ,A B A B ∈ e 1 2, :p p → são funções conhecidas. Note-se que, neste
problema, a coordenada espacial z desempenha o papel do tempo t nas equações das taxas.
Para a resolução deste novo sistema é necessário conhecer as condições iniciais
( )( )
1 0
2 0
00
0p a
zp b
=⎧⎪= → ⎨ =⎪⎩.
Apresenta-se, então, um programa MATLAB que permite resolver este novo problema.
1
02
01
2
1101constantes12escolhidas
1
AaAbB
B
= −==
→==
=
[ ]intervalo de pesquisa 0,10z→ ∈
Programa MATLAB: «diferen.m».
Nota: O programa «diferen» chama, durante a sua execução, a subrotina «derv» através do
programa «ode45».
4 Carlos R. Paiva
Programa «diferen.m»
% DIFEREN Resolução do seguinte sistema de equações diferenciais
%
% d B1 * p1(z)
% ---- p1(z) = ----------------------- * [A1 * p2(z) - 1]
% d z 1 + p1(z) + p2(z)
%
% d B2 * p2(z)
% ---- p2(z) = -------------------- * [A2 * p1(z) - 1]
% d z 1 + p1(z) + p2(z)
%
clear all % Limpa todas as variáveis da memória
close all % Fecha todas as janelas de figuras
A1 = -1; A2 = 1; B1 = 2; B2 = 1; % Constantes do sistema de equações
z=linspace(0,10,1000); % Vector distância onde z = 0 até z = 10 com 1000 pontos
pini=[10 1]; % Valores iniciais das variaveis p1(z) = 10 e p2(z) = 1 no ponto z = 0
OPTIONS=odeset('AbsTol',1e-9,'RelTol',1e-6); % Alteração do valor da tolerância do
% método de RUNGE-KUTTA
% (função ODE45.m)
% Método de RUNGE-KUTTA. 'derv' é a função onde está descrito o sistema diferencial
[z,p]=ode45('derv',z,pini,OPTIONS,A1,A2,B1,B2);
plot(z,p(:,1),'b-',z,p(:,2),'r-'); % Gráfico de p1(z) e p2(z)
title('p1(z) - Azul p2(z) - Vermelho');
xlabel('Distância');
ylabel('Amplitude');
Relatório: Fotónica 5
Função «derv.m»
% DERV Função executada por ode45.m com a descrição do sistema de equações diferenciais
% z - posição onde calcular o valor das derivadas de p1(z) e p2(z)
% p - vector com os valorer p1(z) e p2(z) na posição z
% flag - flags usadas pela função ode45.m
% A1, A2, B1, B2 - constantes do sistema diferencial
%
function dp=derv(z,p,flag,A1,A2,B1,B2)
dp(1)=B1*p(1)/(1+p(1)+p(2))*(A1*p(2)-1); % Equação dp1(z)/dz
dp(2)=B2*p(2)/(1+p(1)+p(2))*(A2*p(1)-1); % Equação dp2(z)/dz
dp=dp'; % Passagem do vector das derivadas de linha para coluna
Figura 1 – Resultado obtido através do programa DIFEREN.
William E. Boyce and Richard C. DiPrima, Elementary Differential Equations and
Boundary Value Problems (Hoboken, New Jersey: Wiley, 8th ed., 2005).
6 Carlos R. Paiva
■
Apresentam-se, de seguida, alguns gráficos correspondentes à simulação numérica das
equações das taxas num laser semicondutor com modulação directa da corrente de injecção.
Relatório: Fotónica 7
Figura 2 – Alguns resultados do Trabalho T1 – alínea (b).
8 Carlos R. Paiva
2. Resolução numérica de equações modais Em vários pontos do programa de Fotónica é necessário resolver numericamente equações
modais. É o caso dos lasers semicondutores, onde é necessário determinar – usando o método
do índice de refracção efectivo – os valores de effn e de n , ver Trabalho T3. É ainda o caso
das fibras ópticas, onde é necessário resolver a respectiva equação modal – seja de forma
rigorosa (modos HE, EH, TE e TM), seja de forma aproximada (modos LP), ver Trabalhos T4
e T5. Em qualquer dos casos a utilização do programa MATLAB «fzero» é perfeitamente
suficiente. O primeiro exemplo que se vai apresentar refere-se ao Trabalho T2 sobre mecânica
quântica. Neste trabalho considera-se um modelo de Kronig-Penney simplificado para
potenciais periódicos. O objectivo é o de entender o aparecimento de bandas de energia
permitidas e proibidas neste tipo de potenciais – uma base indispensável para a física do
estado sólido. Na Fig. 3 apresenta-se um gráfico de (ver Trabalho T2)
( ) ( ) 2 2 2
0 20
,2
k aq a ma
ξ ζ ζ ξ πζ π= ⎡ ⎤
→ = =⎢ ⎥= ⎣ ⎦
EE
E.
Figura 3 – Resultados da segunda parte do Trabalho T2.
Relatório: Fotónica 9
Nas Figs. 4 e 5 apresentam-se os resultados reslativos ao trabalho T3. Nas Figs. 6 e 7
apresentam-se os resultados relativos aos trabalhos T4 e T5 (respectivamente).
Figura 4 – Resultados da primeira fase do Trabalho T3.
George Lindfield and John Penny, Numerical Methods Using MATLAB (London: Ellis
Horwood, 1995).
10 Carlos R. Paiva
Figura 5 – Resultados da segunda fase do Trabalho T3.
Figura 6 – Resultados do Trabalho T4.
Relatório: Fotónica 11
Figura 7 – Resultados do Trabalho T5.
Amnon Yariv and Pochi Yeh, Photonics: Optical Electronics in Modern
Communications (New Yor: Oxford University Press, 6th ed., 2007), pp. 110-155, pp.
797-811.
Bahaa E. A. Saleh and Marvin Carl Teich, Fundamentals of Photonics (Hoboken, New
Jersey: Wiley, 2nd ed., 2007), pp. 289-364.
12 Carlos R. Paiva
3. Propagação de impulsos em fibras ópticas Para a simulação numérica da propagação de impulsos em fibras ópticas há que recorrer à
FFT (fast Fourier transform). Consideremos apenas o caso das fibras monomodais operadas
em regime linear. Despreza-se, aqui, a influência da dispersão de ordem superior (i.e.,
considera-se que é razoável fazer 3 0β = ). Assim, só o efeito da DVG (dispersão da
velocidade de grupo) é tido em consideração. As perdas são desprezadas de forma a poder
analisar, isoladamente, o efeito da DVG sobre a propagação os impulsos. Nestas
circunstâncias e em termos das variáveis normalizadas ζ e τ
220
1
0
D
z zL
t z
βζ
τ
βττ
= =
−=
a equação de propagação dos impulsos é dada pela equação diferencial
( )2
2 2
1 sgn 02
u ui βζ τ∂ ∂
− =∂ ∂
.
Consideremos o caso específico da zona de dispersão anómala em que
( )2
zona de dispersãosgn 1
anómalaβ→ = −
de forma que a equação de propagação dos impulsos se reduz a
2
2
12
u uiζ τ∂ ∂
= −∂ ∂
.
Esta equação é facilmente resolvida no domínio da frequência (normalizada) ξ
( ) ( ) ( ), , ,u u uζ τ ζ ξ ζ τ= ⎡ ⎤⎣ ⎦F
Relatório: Fotónica 13
( ) ( ) ( )
( ) ( ) ( )
, , exp
1, , exp2
u u i d
u u i d
ζ ξ ζ τ ξ τ τ
ζ ξ ζ τ ξ τ ξπ
∞
−∞
∞
−∞
=
= −
∫
∫
uma vez que
( )21 ,2
u i uξ ζ ξζ∂
=∂
( ) ( ) 21, 0, exp2
u u iζ ξ ξ ξ ζ⎛ ⎞∴ = ⎜ ⎟⎝ ⎠
.
Para determinar a forma do sinal, após uma determinada distância de propagação ao longo da
fibra óptica, há que voltar para o domínio do tempo – aplicando agora uma IFFT (a inversa de
uma FFT). No programa que a seguir se apresenta, intitulado «propag», apresenta-se o efeito
da DVG sobre um impulso com uma forma inicial
( ) ( ) ( )0 0, sechu uτ τ τ= = .
Sublinhe-se a necessidade de alterar a ordem interna do vector das frequências devido à forma
como a FFT apresenta a função no domínio da frequência. Este mesmo efeito poderia ser
alcançado através da instrução «fftshift» que, em vez de alterar a ordem interna do vector das
frequências, altera a própria disposição da transformada (colocando o valor correspondente à
frequência nula no centro). Em qualquer caso é necessário garantir que a multiplicação pela
função de transferência da fibra óptica se processa, no domínio da frequência, de forma
correcta.
Amnon Yariv and Pochi Yeh, Photonics: Optical Electronics in Modern
Communications (New Yor: Oxford University Press, 6th ed., 2007), pp. 317-322.
14 Carlos R. Paiva
% PROPAG Propagação de impulsos regida pela equação diferencial
%
% d i d ^ 2
% ---- a(z, t) = --- ---------- a(z,t)
% d z 2 d t ^ 2
%
clear all; % Limpa todas as variáveis da memória
ztotal=5; % Distância total considerada
Nz=500; % Número de passos na distância
z=linspace(0,ztotal,Nz); % Vector distância com 'Nz' posições
t0=20; % Limite da janela temporal
Nt=1024; % Número de amostras temporais
t=linspace(-t0,t0,Nt); % Vector temporal com 'Nt' pontos (Janela temporal)
Ts=t(2)-t(1); % Separação entre amostras
Ws=2*pi/Ts; % Largura total da janela espectral
W=Ws*[0:1:Nt/2-1 -Nt/2:1:-1]/Nt; % Vector de frequências [rad/s]
a=sech(t); % Impulso inicial: a(0,t)=sech(t)
A=fft(a); % Impulso no domínio espectral: A(0,w)
for i=2:Nz % Ciclo para cada posição de z
Az=A.*exp(-1i/2*W.^2*z(i)); % Impulso no domínio espectral: A(z,w)
a(i,:)=ifft(Az); % Impulso no domínio temporal: a(z,t)
end
figure(1); % Gráfico 2D do impulso inicial/final
plot(t,a(1,:),'b',t,abs(a(end,:)),'r');
title('entrada - azul saída - vermelho');
xlabel('Tempo');
ylabel('Amplitude');
axis([-t0 t0 0 1]);
[X,Y]=meshgrid(t,z); % Cria uma grelha rectangular de dimensões (Nt,Nz) usada por 'mesh'
figure(2);
Relatório: Fotónica 15
mesh(X,Y,abs(a)); % Gráfico da evolução do impulso
xlabel('Tempo');
ylabel('Distância');
zlabel('Amplitude');
figure(3);
mesh(X,Y,abs(a));
xlabel('Tempo');
ylabel('Distância');
zlabel('Amplitude');
view(30,45); % view (azimute, elevação);
Figura 8 – Primeira figura obtida através do programa PROPAG.
16 Carlos R. Paiva
Figura 9 – Segunda figura obtida através do programa PROPAG.
Figura 10 – Terceira figura obtida através do programa PROPAG.
Relatório: Fotónica 17
4. Propagação de solitões em fibras ópticas Para estudar a propagação de impulsos numa fibra óptica em regime não linear usa-se o
SSFM (split-step Fourier method). Em termos das variáveis normalizadas ζ e τ , tem-se
então
( ) ( )2
22 2
equação de propagação de impulsos 1 sgnregime não linear 2 2
u ui u u i uβζ τ∂ ∂ Γ
→ − + = −∂ ∂
que se reduz à chamada equação NLS (não linear de Schrödinger) quando ( )2sgn 1β = −
(zona de dispersão anómala) e 0Γ = (sem perdas). Note-se que os efeitos não lineares e
dispersivos de ordem superior são, aqui, desprezados. O SSFM baseia-se na separação entre
os efeitos não lineares (operador N ) e os efeitos dispersivos (operadores Dτ ou Dξ ):
( ) ( )( )
( )
2
2 22
22
1 sgn2 1, sgn
22
D iu D N u D i
N i u
τ
τ ξ
βτ
ζ τ β ξζ
⎧ ∂=⎪ ∂∂ ⎪= + → → = −⎨∂ Γ⎪ = − +⎪⎩
( ) ( ) ( ), , ,u u uζ τ ζ ξ ζ τ= ⎡ ⎤⎣ ⎦F .
Assim, nesta versão do SSFM, considera-se o processo iterativo de passo h que permite ir do
impulso ( ) ( )0 0,u uτ τ= até à saída ( ),Lu ζ τ , com L DL Lζ = .
( ) ( ) ( )SSFM : , , ,u u h wζ τ ζ τ ζ τ→ + =
( ) ( ) ( )( ) ( )( ) ( ) ( )( ) ( )
, exp ,
, FFT ,
, exp ,
, IFFT ,
v h N u
V v
W h D V
w W
ξ
ζ τ ζ τ
ζ ξ ζ τ
ζ ξ ζ ξ
ζ τ ζ ξ
=
= ⎡ ⎤⎣ ⎦
=
= ⎡ ⎤⎣ ⎦
18 Carlos R. Paiva
% SOL_N programa para representar a evolução de um SOLITÃO de ordem N
% Método: SSFM (Split-Step Fourier Method)
%
clear all;
close all;
N=3; % Ordem N do solitão sob consideração
Ntau=2^11; % Número de pontos da janela temporal
tmax=5; % Limites máximo e mínimo da janela temporal
tau=linspace(-tmax,tmax,Ntau); % Definiçao da janela temporal (em 'tau')
zetaMax=pi/2; % Definiçãoo da distância normalizada 'zeta'
h=0.005; % Passo incremental do SSFM
Nstop=round(zetaMax/h); % Número de iterações do SSFM
zeta=linspace(0,zetaMax,Nstop); % Definição da janela espacial (em 'zeta')
Ts=tau(2)-tau(1); % Discretização temporal
Ws=2*pi/Ts; % Discretização na frequência
xi=Ws*[0:1:Ntau/2-1 -Ntau/2:1:-1]/Ntau; % Vector das frequências normalizadas
u=N*sech(tau); % Definição do impulso inicial
beta2=-1; % Sinal de beta 2 (DVG)
gama=0; % Perdas normalizadas
impulso(1,:)=abs(u); % Módulo do impulso
D=i/2*sign(beta2)*xi.^2; % operador D de dispersão (DVG)
for k=2:Nstop
N=-gama/2+i*abs(u).^2;
v=(u).*exp(h*N);
V=fft(v);
U=V.*exp(h*D);
u=ifft(U);
impulso(k,:)=abs(u);
end
Relatório: Fotónica 19
[XN,YN]=meshgrid(tau,zeta);
mesh(XN,YN,impulso);
axis([-tmax tmax 0 zetaMax 0 max(max(impulso))]);
xlabel('\tau = (t - \beta_1 z) / \tau_0');
ylabel('\zeta = z / L_D');
zlabel('|u| (\zeta, \tau)');
title('Solitão de ordem N = 3')
Figura 11 – Resultado obtido através do programa SOL_N.
Govind P. Agrawal, Nonlinear Fiber Optics (San Diego, California: Academic Press,
4th ed., 2007), pp. 41-45.
Recommended