PROPRIEDADES ELECTROMAGNÉTICAS DOS MATERIAIS
1ª Série
A - SIMULAÇÃO NUMÉRICA DO TRANSPORTE DE CARGA
IST, Fevereiro 2006
TRANSP1
SIMULAÇÃO NUMÉRICA DO TRANSPORTE DE CARGA
I - INTRODUÇÃO
A simulação numérica não substitui a realidade física. Todavia, quando os modelos utilizados
são convenientemente validados por confronto com resultados experimentais, a simulação
pode contribuir para a compreensão dos fenómenos físicos.
A corrente eléctrica é o resultado do movimento de partículas com carga, genericamente
designado transporte de carga eléctrica. Esse transporte assume essencialmente duas formas:
i)
•
Por acção de uma força exterior;
•
Por difusão.
Um exemplo da acção de uma força exterior é a aplicação de uma tensão aos terminais
de uma resistência (metálica, de semicondutor, etc.). Devido a essa acção exterior
aparece no interior do material um campo eléctrico e, portanto, uma força eléctrica. Esse
G
transporte designa-se por condução e a relação entre a densidade de corrente J e o
G
campo eléctrico E é dada pela lei de Ohm na forma local:
G
G
J = σE
(1)
onde σ é a condutividade. Subjacente ao resultado anterior está que a densidade de
G
corrente associada a partículas com velocidade v é dada por:
G
G
J = ρv
(2)
onde ρ é a densidade volumétrica de carga. A velocidade referida é de facto o valor
médio da velocidade das partículas. Esse valor médio é zero na ausência de força
exterior, isto é, em equilíbrio termodinâmico, embora em valor instantâneo a velocidade
das partículas não seja nula. À temperatura ambiente, na ausência de campo, a
velocidade das partículas com massa eficaz igual à massa do electrão em repouso é da
ordem de 105m/s. Contudo, segundo qualquer direcção o número de partículas que se
deslocam num sentido é igual ao número das que se deslocam em sentido contrário. É
desta “desorganização” total que resulta uma velocidade média nula. A hipótese
ergódica é válida na situação estacionária. Diz-nos que a média tomada sobre a
TRANSP2
população de portadores de carga num dado instante conduz ao mesmo resultado que a
média tomada sobre um portador ao longo do tempo. Na ausência de campo, o
movimento desordenado de uma partícula corresponde a percursos rectos, com
velocidade constante, que são interrompidos por “colisões” com imperfeições da rede
cristalina. Estas fazem mudar a direcção da partícula de forma aleatória1 , isto é, “à
sorte”, ao acaso, significando que a direcção após a colisão pode ser qualquer. O campo
eléctrico vem transformar as trajectórias entre colisões, que eram segmentos de recta, em
arcos de parábola. O valor médio da velocidade deixa de ser zero, passando a ser para
campos fracos dado por
G
G
v = ±μ E
(3)
onde μ é a mobilidade dada por:
μ=
Q τL
mef
(4)
onde Q é a carga da partícula, mef é a massa eficaz da partícula e τ L é o tempo que em
média as partículas estão sem colidir com as imperfeições da rede cristalina e que se
chama tempo de livre percurso médio. A presença do módulo na carga deve-se a que as
mobilidades são sempre positivas. Em (3) o sinal depende da carga do portador, devendo
usar-se o sinal positivo (negativo) para partículas de carga positiva (negativa).
Subjacente à noção de livre percurso médio está a forma como as partículas colidem com
as imperfeições ao longo do tempo. Se considerarmos que em t=0 há N0 partículas que
vão começar a colidir, no instante t não terá ainda colidido um número de partículas N
dado por
N (t ) = N 0 exp( −t / τ L )
(5)
Por diferenciação do resultado anterior o número de partículas Δn que no intervalo de
tempo Δt colide é dado por
Δn = − n
Δt
τL
(6)
onde o sinal menos mostra que a população que não colidiu decresce.
A contribuição para a condutividade devida a partículas com concentração N na unidade
1
Do latim “alea”, sorte.
TRANSP3
de volume é Q N μ . Nos metais e ligas metálicas temos electrões. Nos semicondutores
somam-se as contribuições de electrões e de buracos, sendo a condutividade dada por:
(
σ = Q n μn + p μ p
)
(7)
para concentrações volumétricas de electrões e de buracos dadas por n e p,
respectivamente. Se a força exterior além de eléctrica fosse também magnética tínhamos
o chamado efeito de Hall, que pode ser forte nos semicondutores mas é fraco nos metais.
ii) Outra forma de transporte de carga é a difusão, a que corresponde a corrente de difusão.
É o processo pelo qual se dá o “espalhamento” das partículas das concentrações mais
altas para as mais baixas. É o fenómeno que justifica a uniformização da temperatura
num corpo que seja aquecido numa extremidade ou o espalhamento de uma gota de tinta
que se deite na água.
Para haver difusão é assim preciso que haja uma concentração não uniforme, isto é um
gradiente da concentração das partículas e agitação térmica. Esta agitação térmica
garante precisamente a desordem que, num modelo simplificado a uma dimensão de
espaço, faz com que a partícula ao colidir com uma imperfeição tanto possa ir para a
frente como para trás, com igual probabilidade. O movimento que corresponde a este
comportamento é chamado movimento Browniano. A equação que rege a difusão é a
chamada equação de difusão ou do calor que, para uma dimensão de espaço físico, é
dada por:
∂N
∂2 N
=D 2
∂t
∂x
(8)
onde t é o tempo, x a coordenada de espaço e D o coeficiente de difusão.
A densidade de corrente eléctrica de difusão é dada por:
G
J = −Q D grad N
(9)
A difusão só é possível nos semicondutores, visto que nos metais não é possível ter uma
concentração de portadores não uniforme. Em (9) o sinal “−“ resulta do facto de a
G G
corrente de partículas por difusão C = J Q se dar sempre no sentido decrescente da
concentração de portadores, ou seja, no sentido contrário ao gradiente de partículas.
TRANSP4
II - APLICAÇÕES
Atendendo a que o comportamento aleatório é essencial quer para a condução quer para a
difusão, os programas de simulação numérica que vão ser utilizados no presente trabalho
utilizam um gerador de números pseudo-aleatórios. A designação resulta do facto de a
sequência de números gerados ser definida deterministicamente a partir do momento em que
se gera de forma aleatória o primeiro número da série. A gama de números gerados deve
garantir a função de distribuição que se pretende simular. As decisões sobre a possibilidade
de “colisão” e as alterações no movimento da partícula em caso de “colisão” serão tomadas
com base nos números gerados pseudo-aletaoriamente. Os resultados são obtidos
“observando” sucessivas trajectórias das partículas simuladas (ou só uma ao longo do tempo
se a situação for estacionária), alternando-se deriva com colisões pseudo-aleatórias e tomando
a média sobre elas. Trata-se assim de uma média sobre sucessivas “experiências simuladas”.
Estes métodos são também designados por Métodos de Monte Carlo, e dado que o
comportamento não é determinístico, há flutuações nos resultados. Realizar mais
“experiências”, isto é, tomar médias sobre populações maiores, reduz as flutuações. É
contudo possível mostrar pelo teorema do limite central que, atendendo ao carácter gaussiano
da distribuição, as flutuações decrescem com a raiz quadrada do número de “experiências”
(variância inversamente proporcional ao número de “experiências”). Assim, se fizermos 100
vezes mais experiências, o que multiplica por 100 o tempo de cálculo, o desvio quadrático
médio será só 10 vezes menor.
Note-se que por razões de tempo de cálculo as aproximações utilizadas no programa para
simulação da condução e da difusão são diferentes, uma vez que para a difusão o intervalo de
tempo de análise Δt coincide com o tempo de livre percurso médio τL.
II.A – Transporte por acção de uma força exterior
Programa: CONDF. PAS;
Linguagem: PASCAL7
Neste programa o comportamento de uma partícula é analisado em sucessivos intervalos de
tempo Δt. Durante um intervalo de tempo Δt, admitindo que ao campo electrico corresponde
uma força aceleradora segundo x, as velocidades variam de acordo com
v x = v0 x +
qE
m*
Δt
v y = v0 y
(10)
TRANSP5
e as coordenadas de espaço como
Δx =
qE 2
Δt + v0 x Δt
2m*
Δy = v0 y Δt
(11)
A probabilidade deste intervalo terminar com uma colisão é
Δt
τL
(12)
sendo o gerador de números pseudo aleatórios que decide se há ou não colisão. Caso haja é
também ele que decide a nova direcção da velocidade. Todas as direcções são equiprováveis
e o módulo da velocidade mantém-se. Execute o programa segundo as instruções e crie os
ficheiros. Usando o MATLAB crie as variáveis a partir dos ficheiros e faça os gráficos (ver
Apêndice).
II.A.1 – Agitação térmica ( E = 0 , força nula)
1. Com campo nulo, em que a velocidade da partícula resulta apenas da agitação
térmica, iniciar o gerador de números pseudo-aleatórios. Para isso digite
97,98…102.
2. Tome nota do valor médio da velocidade segundo x. Compare os resultados com
o valor que obteria de (3). Tome como gerador o que der menor módulo da
velocidade segundo x.
3. Obtenha as trajectórias (x, y) e as evoluções temporais das componentes das
velocidades vx (t ) e v y (t ) .
II.A.2 – Condução ( E ≠ 0 )
1. Fixado o gerador, estude as situações correspondentes a E = 105 , 5 × 105 e
106 V/m. Calcule o valor da mobilidade dada pela relação (3) e verifique que,
quando o campo cresce, a mobilidade tende para o valor teórico dado por (4).
Porquê?
2. Observe as trajectórias (x, y) e as evoluções temporais das componentes das
velocidades vx (t ) e v y (t ) para E = 106 V / m . Comente estes resultados e
compare-os com os obtidos para E = 0 .
3. Repetir a experiência 1., mas com um gerador a que corresponda uma velocidade
de agitação térmica de sinal contrário ao que tomou em 1.
TRANSP6
II.B – Transporte por acção difusão
II.B.1 – Difusão a uma dimensão partindo de um Dirac em x = 0
Programa: DIF1D. PAS;
Linguagem: PASCAL7.
Considere uma distribuição inicial de N 0 partículas correspondente a um impulso de Dirac
centrado em x = 0 . A solução de (8) para esta situação está associada a uma evolução no
tempo e no espaço dada por [1]:
−
N0
N ( x, t ) =
2πσ
2
e
x2
2σ2
(13)
A equação (13) corresponde a uma função densidade gaussiana de valor médio nulo, sendo a
variância dada por:
σ2 = ( x − x ) = x 2 − x 2 = 2 Dt
2
(14)
Vemos assim que a variância cresce linearmente no tempo.
A simulação é feita a uma dimensão de espaço, dividindo o espaço em análise x em intervalos
Δx e o tempo de análise t em intervalos Δt.
Parte-se de uma situação inicial em que todas as partículas estão em x=0.
Admite-se simplificadamente que as partículas têm todas o mesmo módulo da velocidade
G
v0 ,mas que esta pode ser segundo +x ou –x com igual probabilidade. Em cada intervalo Δt as
G
partículas deslocam-se de v0 Δt. Este intervalo de tempo termina com uma “colisão” em que
com base num gerador de números pseudo-aleatório se decide com igual probabilidade se o
sentido da velocidade se mantem ou troca.
Ao fim de cada Δt são contadas as partículas em cada intervalo Δx e obtem-se nesse instante a
distribuição ao longo de x. Esta distribuição tem um perfil aproximadamente gaussiano (as
condições do teorema do limite central verificam-se aproximadamente) que vai evoluindo de
Δt em Δt com progressivo “alargamento” da distribuição.
Para o movimento Browniano [2] o coeficiente de difusão correspondente é dado por
D = v02Δt 2
(15)
TRANSP7
Os resultados só fazem sentido em “média”, isto é, considerando muitas partículas e tomando
o valor médio do seu comportamento. Quanto mais são as partículas utilizadas na simulação,
menores as flutuações e melhores os resultados, mas também maiores são a memória
necessária e o tempo de cálculo.
Por limitações da linguagem, o número de variáveis é limitado, sendo também limitado o
número de partículas que se pode simular. Para manter as flutuações dentro de valores
aceitáveis, e apesar da situação não ser estacionária, estão a tomar-se também médias no
tempo. Desde que a evolução seja lenta, os valores são representativos do comportamento no
instante médio desse intervalo.
Procedimento:
Execute o programa seguindo as instruções que este lhe der. Irá obter ficheiros com duas
distribuições relativas a dois “instantes” diferentes que poderão ser visualizadas e impressas
utilizando o programa MATLAB. No programa crie as variáveis a partir dos ficheiros e faça
os gráficos (ver Apêndice). A partir das distribuições obtenha a evolução da variância no
tempo e calcule o coeficiente de difusão usando (13) e (14). Compare com os resultados que
o programa calcula directamente pelas flutuações e com os que pode obter de (15), admitindo
que v0 = 105 m / s e Δt = 1 ps .
II.B.2 – Difusão a duas dimensões partindo de um impulso de Dirac
Programa: DIF2D. PAS;
Linguagem: PASCAL7.
É descrito pela equação:
⎛ ∂2 N ∂2 N ⎞
∂N
= D⎜ 2 + 2 ⎟
⎜ ∂x
∂t
∂y ⎟⎠
⎝
(16)
Do ponto de vista do programa, o tratamento é inteiramente análogo ao feito a uma dimensão.
Contudo, agora há duas direcções, x e y tratadas cada uma delas através de intervalos Δx e Δy,
utilizando-se o gerador de números pseudo-aleatórios, independentemente para cada uma das
direcções, para decidir do sentido após as colisões.
TRANSP8
Procedimento:
Siga as instruções do programa DIF2D.pas para obter a evolução temporal da distribuição de
partículas a duas dimensões de espaço físico. Trace os gráficos de n(x,y) e y(x) para o
primeiro e último instantes de tempo.
II.B.3 – Difusão com recombinação
Programa: DIFREC1F. PAS;
Linguagem: PASCAL7.
A equação diferencial não ordinária que corresponde ao problema a uma dimensão de espaço
é
N − N0
∂N
∂2 N
=−
+D 2
∂t
τv
∂x
(17)
onde τv é o tempo de vida médio (representa o tempo que em média uma partícula “vive”
quando não existe o termo dependente de x), N 0 é o valor da concentração de partículas na
situação estacionária e a segunda parcela do segundo membro está ligada com o mecanismo
de transporte por difusão. Na situação estacionária a equação é
N − N0
d 2N
+D 2
0=−
τv
dx
(18)
Para cristais semi-infinitos e para uma concentração fixa em x = 0 , N ( x = 0) = Nin , obtém-se
o andamento:
(
N ( x ) = N 0 + ( Nin − N 0 ) exp − x Ldif
)
(19)
onde Ldif é o comprimento de difusão dado por
Ldif = D τv
(20)
Este programa difere do anterior porque a concentração de partículas numa fronteira dada é
imposta com um valor constante no tempo e porque as partículas desaparecem por
recombinação. Vamos supor que a outra fronteira está a uma distância muito maior do que o
comprimento de difusão e tem uma concentração de partículas nula.
TRANSP9
Um exemplo deste comportamento é a difusão de portadores de minoria num semicondutor,
como é o caso dos electrões do lado p de um díodo (ou de buracos do lado n), a partir da zona
de transição em direcção ao contacto.
A análise no tempo e no espaço é feita como no programa anterior em termos de intervalos Δt
e Δx dados.
Seja τv o tempo de vida, ou seja o tempo que as partículas em média demoram a
“desaparecer” por recombinação. No intervalo de tempo Δt a probabilidade de cada uma das
partículas se recombinar é Δt/τv . É também o gerador pseudo-aleatório que ao fim de cada
intervalo Δt decide se as partículas se recombinam ou não respeitando a probabilidade
anterior.
No fim de cada intervalo de tempo a concentração pretendida na fronteira (a primeira “fatia”
de Δx) é também reposta, retirando partículas aleatoriamente se estão a mais, ou
introduzindo-as, com velocidades com igual probabilidade segundo +x e –x utilizando
também o gerador pseudo-aleatório, se estão a menos.
Se a situação for estacionária, além de valores médios sobre a população de partículas podem
tomar-se valores médios ao longo de vários Δt. Verifica-se a hipótese ergódica. No programa
deixa-se passar um tempo de simulação muito maior do que o tempo de vida médio, para que
a situação estacionária seja atingida, tomando-se depois médias sobre as partículas ao longo
do tempo.
No caso presente o valor de N0 em (18) é nulo. Sendo xi a coordenada do ponto médio do
intervalo i de largura Δx e Ni o número de partículas nesse intervalo, temos para a distribuição
espacial de partículas
N i = N 1e
−
( xi − x1 )
Ldif
(21)
onde o comprimento de difusão é dado por
Ldif = v0 Δtτ v
(22)
As velocidades das partículas dependem da temperatura e da massa eficaz, portanto do
material. São correntes massas eficazes entre 10 vezes maiores e 10 vezes menores do que a
do electrão em repouso. Considerando um electrão em repouso e uma distribuição à Maxwell
Boltzmann teríamos, a três dimensões de espaço, v0 ≅ 105 m / s .
TRANSP10
Os tempos de percurso livre de colisões, isto é, entre colisões, variam muito. Dependem do
material, da temperatura, da energia das partículas, da concentração de partículas, etc.. A
300 K podem tomar muitos valores, sendo correntes valores de 10-13 a 10-15s.
Os tempos de vida são geralmente maiores do que os de livre percurso. Há no entanto
excepções (é o caso por exemplo no “spectral hole burning” 2 dos laser). Tal como os tempos
de percurso livre, os tempos de vida dependem do material, da temperatura, da energia das
partículas, da concentração de partículas, etc.. A 300 K podem tomar muitos valores, sendo
correntes valores de 10-10 a 10-13s.
Ao escolher valores para o programa há várias limitações a ter em conta. Por um lado deverá
ser Δx << Ldif ou seja Δx << v 0 Δtτ v para que a distribuição espacial das partículas fique
bem definida.
Por outro lado deverá ser Δx >> v0 Δt , para que as partículas não atravessem um ou mais
intervalos Δx sem serem contabilizadas.
As relações anteriores tornam desejável que
Δt
τv
<< 1
(23)
Baixar Δt faz subir o tempo de cálculo para a simulação de um mesmo tempo físico t. O
mesmo acontece com aumentar τv , já que o tempo de análise lhe deve ser muito superior para
que seja atingida a situação estacionária.
O aumentar τv faz ainda aumentar o número de partículas envolvidas, já que a distribuição se
“espalha” mais no espaço, por aumento do comprimento de difusão, e aumentando por isso a
memória e tempo de cálculo necessários. O aumentar v0 tem um efeito semelhante mas mais
forte, de acordo com a definição de comprimento de difusão.
Procedimento:
Obtenha os ficheiros de valores correspondentes às opções 1 a 4. O programa cria para cada
opção dois ficheiros correspondentes a instantes de tempo diferentes, quer para a posição x,
quer para a densidade de partículas n. Nestas condições:
a) Verifique que em cada opção é atingida a situação estacionária, comparando as curvas n(x)
correspondentes a instantes diferentes.
2
Designação anglo-saxónica para o desaparecimento da emissão em determinadas frequências com
o aumento da potência emitida.
TRANSP11
b) Trace no mesmo gráfico as curvas n(x) correspondentes ao instante de tempo final das
opções anteriores verificando:
b.1) A influência de v0, comparando as opções 1 e 2.
b.2) A influência do tempo de vida das partículas, comparando as opções 1 e 3.
b.3) A influência do tempo de colisão das partículas, comparando as opções 1 e 4.
Nota: Ficheiros criados para a opção i e instante de tempo j, xij e nij.
Programa: DIFREC2F. PAS;
Linguagem: PASCAL7.
O programa é inteiramente análogo ao programa anterior, DIFREC1F. A diferença reside em
que é possível fixar as concentrações em duas fronteiras e variar o afastamento entre elas. A
versão DIFREC1F, que é um caso particular deste programa, existe por razões de eficiência,
visto que só com uma fronteira com concentração imposta é possível reduzir o tempo de
cálculo significativamente.
Procedimento:
Corra as opções que o programa oferece, correspondentes a diferentes afastamentos entre as
fronteiras.
Nota: A opção i cria os ficheiros xi e ni.
III - BIBLIOGRAFIA
[1] A Sommerfeld, “Partial Differential Equations in Physics”, Academic Press, 1949.
[2] Rozanov, “Processus Aléatoires”, Cap. II, pag. 85 e seguintes, MIR, 1975.
TRANSP12
APÊNDICE – INSTRUÇÕES MATLAB
Exemplo: Ficheiros x (posições) e n (concentrações)
load x cria o vector x
;
load n cria o vector n
plot (x, n) cria o gráfico n(x)
plot (x1, n1, x2, n2) cria os gráficos simultâneos n1(x1) e n2(x2)
É possível especificar a cor, o tipo de linha e os símbolos quando se desenha o gráfico com o
comando plot.
plot (x, n, ’cor_estilo_símbolo’)
Cores: ‘c’ – ciânico (azul)
‘m’ – magenta
‘y’ – amarelo
‘r’ – vermelho
‘g’ – verde
‘b’ – azul
‘w’ – branco
‘k’ – preto
Estilos: ‘-‘ – cheio
‘- -‘ – tracejado
‘:’ – ponteado
‘- .’ – traço ponto
omissão de estilo não considera a linha
Símbolos: ‘+’
‘o’
‘*’
‘×’
‘s’ – quadrado
‘d’ – diamante
‘^’ – triângulo
Exemplo: plot (x, n, ‘ks’) – traça n(x) com quadrados pretos em cada ponto mas não liga
os pontos por linhas.
plot (x, n, ‘r:+’) – traça n(x) com uma linha vermelha a ponteado e põe símbolos
+ em cada ponto.
TRANSP13
IV – QUESTIONÁRIO
ESTE QUESTIONÁRIO DESTINA-SE A QUE OS ALUNOS FAÇAM UM BALANÇO DOS SEUS
CONHECIMENTOS SOBRE ESTA MATÉRIA. SE TIVER DÚVIDAS TIRE-AS NUM DOS
HORÁRIOS DE DÚVIDAS. NÃO INCLUA AS RESPOSTAS NO RELATÓRIO.
1. Considere duas resistências de material intrínseco, uma de germânio e outra de silício.
Para cada uma das resistências calcule a variação da temperatura necessária para que o
valor das resistências passe a n vezes o valor que têm a 300 K.
Si (300 K) – WG = 1,12 eV ; Ge (300 K) – WG = 0,66 eV ; μ n , μ p ∝ T −3 / 2
2. Considere duas resistências de semicondutor de igual forma e dimensões, sendo uma de
Germânio e outra de Silício. Para T = 300 K elas tomam o mesmo valor, sendo uma delas
de semicondutor intrínseco.
a) Diga, justificadamente, qual das duas é de material intrínseco e qual a concentração
de impurezas de substituição na outra, que se admite ser do tipo n.
Ge (300 K):
ni = 2, 4 ×1019 m−3 ;
Si (300 K):
ni = 1,5 ×1016 m−3
μ n = 0, 4 m 2V −1s −1 ;
μ n = 0,15 m 2V −1s −1
μ p = 0, 2 m2V −1s −1 ;
μ p = 0, 05 m2V −1s −1
WG = 0, 66 eV
. .
WG = 1,12 e.V .
b) Considere a resistência do material extrínseco iluminada uniformemente de tal modo
que Gfe = m Gtérmico. Suponha que em t = 0 se interrompe abruptamente a iluminação.
Diga ao fim de quanto tempo a concentração de buracos difere menos de b% do seu
valor de equilíbrio. Tempo de vida médio dos portadores de minoria é 10 μs
DADOS:
n
1/2
1/3
1/4
m
4
9
2
b (%)
10
20
5