( PDF ) Estudo do alcance de elétrons com energias entre 110 eV e 50 keV

Download PDF

ads:

“Estudo do alcance de elétrons com energias

entre 110 eV e 50 keV”

Geórgia Santos Joana

Disserta

ão apresentada como parte dos requisitos

para obten

ão do

rau de Mestre em Ciência e

Tecnologia das Radiações, Minerais e Materiais

2005

ads:

Livros Grátis

http://www.livrosgratis.com.br

Milhares de livros grátis para download.

Comissão Nacional de Energia Nuclear

CENTRO DE DESENVOLVIMENTO DA TECNOLOGIA NUCLEAR

Programa de Pós-Graduação em Ciência e Tecnologia das

Radiações, Minerais e Materiais

ESTUDO DO ALCANCE DE ELÉTRONS COM

ENERGIAS ENTRE 110 eV E 50 keV

Geórgia Santos Joana

Dissertação apresentada ao Curso de Pós-Graduação em

Ciência e Tecnologia das Radiações, Minerais e

Materiais, como requisito parcial à obtenção do Grau de

Mestre.

Área de concentração: Aplicação de Técnicas Nucleares

Orientadora: Suely Epsztein Grynberg

Belo Horizonte

2005

ads:

DEDICATÓRIA

Aos meus pais, Jorge e Graça, aos meus

irmãos, Graciela e Giordano e a todos que me

incentivaram e ajudaram de alguma forma a

superar mais esta etapa da minha vida.

iii

AGRADECIMENTO

À Deus, por mais esta graça alcançada. Aos

meus pais, meus irmãos e ao Felipe, por tudo!

(vocês sabem...). Aos meus orientadores (todos

eles). Às pessoas que acreditaram em mim e me

deram a maior força. Ao CDTN e à CNEN.

RESUMO

Existe um grande número de estudos que examina a dependência entre energia e alcance

para elétrons com energias desde alguns keV’s até centenas de MeV’s. Porém, a descrição

quantitativa do transporte de elétrons de baixa energia é um problema complexo,

especialmente devido à não disponibilidade de dados confiáveis para as seções de choque

dos diferentes mecanismos por meio dos quais elétrons interagem com a matéria. Uma vez

que os modelos analíticos de transporte de elétrons ou são incompletos ou apresentam

dificuldades até o momento intransponíveis, os códigos de Monte Carlo têm sido a

alternativa adotada nestas abordagens e são hoje extensivamente utilizados.

É sabido que na interação de radiações ionizantes com o tecido vivo podem ocorrer danos

severos em estruturas nos níveis celular e sub-celular. Uma vez que a transferência de

energia de radiações ionizantes para o meio ocorre predominantemente através da

interação dos elétrons secundários produzidos por esta com o material, existe um interesse

especial na determinação da relação entre energia inicial e alcance para elétrons;

especialmente os de baixa energia.

Neste trabalho é estudada a relação entre energia inicial e alcance em meio biológico para

elétrons com energias entre 110 eV e 50 keV. Para tal foi utilizado o código de Monte

Carlo PENELOPE com o qual foram realizadas simulações do transporte de elétrons de

baixa energia (110 eV a 50 keV) – e dos fótons secundários gerados por estes –

provenientes de uma fonte monoenergética isotrópica localizada no interior de um meio

tecido equivalente homogêneo. A partir daí, as grandezas de interesse foram calculadas e

analisadas.

ABSTRACT

There are a great number of studies that examines the dependence between energy and

range for electrons with energies from some keV's to hundreds of MeV's. However, the

quantitative description of low-energy electron transport is a complex problem, especially

due to the fact that there are not reliable data available for cross sections of the different

possible types of electron interactions with matter. Once the analytic models of electron

transport are incomplete or present difficulties to the present insoluble Monte Carlo codes

have been the alternative adopted in these approaches and are used extensively nowadays.

It is known that the interaction of ionizing radiation with tissue may lead to severe

damages in cellular and sub-cellular levels structures. Thus, there is a special interest in the

determination of the relationship between initial energy and range for electrons despite the

fact that the ionizing radiation energy transfer to the medium occurs predominantly

through the interaction of secondary electrons produced in the material - mainly low-

energy electrons.

In this work it is studied the relationship between initial energy and range in biological

medium for electrons with energies between 110 eV and 50 keV. Monte Carlo PENELOPE

code was used. The simulation of transport of low-energy electrons (110 eV to 50 keV)

originated from a monoenergetic isotropic source located inside a homogeneous tissue-

equivalent medium was accomplished, and then the quantities of interest were calculated

and analyzed.

SUMÁRIO

INTRODUÇÃO................................................................................................................................................8

REVISÃO DE LITERATURA.....................................................................................................................11

1 O MÉTODO MONTE CARLO..........................................................................................................11

1.1 ELEMENTOS DA TEORIA DE PROBABILIDADE ......................................................................................12

1.2 MÉTODOS DE SIMULAÇÃO DE VARIÁVEIS ALEATÓRIAS ......................................................................16

1.2.1 Gerador de números aleatórios...............................................................................................16

1.2.2 Método da transformada inversa.............................................................................................17

1.2.3 Método de Walker....................................................................................................................19

1.2.4 Métodos de rejeição.................................................................................................................21

1.2.5 Métodos de composição. Variáveis bidimensionais.................................................................22

1.3 INTEGRAÇÃO MONTE CARLO .............................................................................................................23

2 INTERAÇÃO DA RADIAÇÃO COM A MATÉRIA.......................................................................27

2.1 INTERAÇÃO DE FÓTONS COM A MATÉRIA............................................................................................27

2.1.1 Efeito fotoelétrico.....................................................................................................................29

2.1.2 Espalhamento Compton...........................................................................................................30

2.2 INTERAÇÃO DE PARTÍCULAS CARREGADAS COM A MATÉRIA ..............................................................31

2.2.1 Poder de frenamento................................................................................................................33

2.2.2 Alcance.............. ......................................................................................................................40

3 TÉCNICAS DE MONTE CARLO PARA O TRANSPORTE DE ELÉTRONS E FÓTONS.......44

3.1 MODELOS CLASSE I E CLASSE II.........................................................................................................46

3.2 SIMULAÇÃO DO TRANSPORTE DE RADIAÇÃO ......................................................................................48

3.2.1 Modelo de espalhamento e funções distribuição de probabilidade.........................................48

3.2.2 Geração das trajetórias aleatórias..........................................................................................51

3.2.3 Seleção de um processo físico..................................................................................................54

3.2.4 Efeitos do tamanho do passo da partícula...............................................................................55

3.2.5 Transporte de partículas como um processo Markoviano.......................................................58

3.2.6 Transporte de elétrons.............................................................................................................59

3.2.7 Transporte de fótons................................................................................................................61

3.3 MÉDIAS ESTATÍSTICAS E INCERTEZAS ................................................................................................62

3.4 EXATIDÃO DO CÓDIGO .......................................................................................................................66

3.4.1 Dados de seção de choque.......................................................................................................66

3.4.2 Erros Sistemáticos ...................................................................................................................67

3.4.3 Incertezas estatísticas ..............................................................................................................68

3.4.4 Eficiência Computacional........................................................................................................69

3.5 REDUÇÃO DE VARIÂNCIA...................................................................................................................70

3.5.1 Interação Forçada...................................................................................................................70

vii

3.5.2 Cisão (splitting) e roleta Russa................................................................................................72

3.5.3 Rejeição de alcance .................................................................................................................73

3.5.4 Outros métodos ......................................................................................................................74

METODOLOGIA..........................................................................................................................................75

4 O CÓDIGO PENELOPE.....................................................................................................................77

4.1 INTERAÇÕES DE FÓTONS.....................................................................................................................77

4.1.1 Espalhamento coerente (Rayleigh)..........................................................................................77

4.1.2 Efeito fotoelétrico.....................................................................................................................77

4.1.3 Espalhamento incoerente (Compton).......................................................................................78

4.1.4 Produção de pares...................................................................................................................78

4.2 INTERAÇÕES DE ELÉTRONS E PÓSITRONS............................................................................................79

4.2.1 Colisão elástica ......................................................................................................................79

4.2.2 Colisão inelástica.....................................................................................................................79

4.2.3 Emissão de bremsstrahlung.....................................................................................................79

4.2.4 Aniquilação de pósitrons .........................................................................................................80

4.3 RELAXAÇÃO ATÔMICA .......................................................................................................................80

4.4 MECÂNICA DO TRANSPORTE DE ELÉTRONS E PÓSITRONS....................................................................81

4.5 SIMULAÇÃO MISTA (CLASSE II) DO ESPALHAMENTO ELÁSTICO ..........................................................82

4.5.1 Espalhamento múltiplo. Método da dobradiça aleatória (“random hinge”)..........................82

4.5.2 Perdas de energia suaves. “Continuous-Slowing-Down-Approximation”..............................82

4.5.3 Combinação de espalhamento e perda de energia ..................................................................83

4.6 GERAÇÃO DAS HISTÓRIAS ..................................................................................................................85

4.7 SIMULAÇÃO DO CHUVEIRO.................................................................................................................89

4.8 ESTABILIDADE DO ALGORITMO DE SIMULAÇÃO.................................................................................93

5 O CÁLCULO DOS ALCANCES........................................................................................................95

RESULTADOS............................................................................................................................................100

CONCLUSÕES............................................................................................................................................111

APÊNDICE A – Composição do tecido mole (ICRP)...............................................................................113

APÊNDICE B – Termos, grandezas e unidades relacionados às radiações ionizantes..........................114

APÊNDICE C – Cinemática de colisão......................................................................................................118

INTRODUÇÃO

Em física, o termo radiação refere-se usualmente a partículas e campos que se propagam

(transferindo energia) no espaço (preenchido ou não por matéria). As radiações podem ser

ionizantes ou não ionizantes. Radiações ionizantes são usualmente definidas como aquelas

que podem ionizar a matéria diretamente ou através da ação de alguma radiação

secundária. O processo de ionização ocorre quando uma quantidade suficiente de energia é

transferida para o material arrancando elétrons dos átomos deste. É importante salientar

que, apesar de serem definidas como “ionizantes”, este não é o único – e nem

necessariamente o mais importante – processo através do qual a radiação pode transferir

energia para um material. Um segundo fenômeno importante é a excitação. Neste processo

o átomo é elevado a um estado de maior energia. No processo de desexcitação ocorre a

emissão de raios-X ou de elétron(s) de camadas mais exteriores àquela onde ocorreu a

excitação (elétrons Auger). Ionizações e excitações podem acarretar sérias conseqüências

físicas, químicas e/ou biológicas (quando cabível). Sistemas biológicos, em especial, são

bastante sensíveis a estes eventos. A transferência de energia de radiações ionizantes para o

meio ocorre predominantemente através da interação dos elétrons secundários com o

material (ATTIX, 1968). A interação destes elétrons com o tecido vivo pode provocar

danos em estruturas celulares dependendo da distribuição local da densidade de energia

absorvida. Para se produzir danos no DNA, por exemplo, a deposição de energia deve ser

extremamente localizada, isto é, a distribuição da energia deve se adequar a uma região

cujo volume é da ordem do diâmetro da hélice desta molécula. Neste contexto, deve-se

fazer referência aos elétrons de baixa energia, como os elétrons Auger, cujos alcances são

tais que se torna possível obter valores altos para a densidade de energia em pequenos

volumes.

Encontra-se na literatura vários exemplos que mostram que elétrons de baixa energia

provenientes de emissores localizados internamente em meios biológicos podem produzir

danos relevantes às entidades biológicas irradiadas desde que a densidade de energia

localmente absorvida seja suficientemente grande (MALAMUT, 1987). Com base nisto,

cogita-se fazer uso de emissores de elétrons Auger, tais como

Br,

201

Tl,

Se,

dentre outros, na radioimunoterapia (BLOOMER, 1978 citado por MALAMUT, 1987;

MAUSNER, 1993). A idéia central desta aplicação é utilizar os elétrons gerados através do

efeito Auger na irradiação de células neoplásicas. Devido à baixa energia destes elétrons

(em geral menor do que 50 keV) devem ocorrer depósitos de energia em regiões

extremamente localizadas. Tais depósitos podem levar a danos severos nas estruturas

moleculares em torno dos locais de decaimento. Se estas estruturas forem vitais para a

viabilidade celular, um efeito letal pode ocorrer (CHARLTON, 1978; TISLJAR-

LENTRILIS, 1978; KASSIS, 1980 citado por MALAMUT, 1987). O limite de 50 keV

para análise detalhada de distribuição de dose nos níveis celular e sub-celular é

corroborado por outros autores ao considerarem que, com energias maiores que 50 keV, os

elétrons, e desta forma a energia, se difundem para longe da fonte (HUMM, 1993).

Portanto, o estudo dos efeitos da radiação no tecido vivo fundamenta-se na análise dos

elétrons de baixa energia. A avaliação do “alcance” destes elétrons se conFIGURA como

uma ferramenta útil na determinação da extensão sobre a qual estes elétrons, colocados em

movimento em torno da trajetória da radiação, induzem danos diretos ao DNA de células

vivas (MEESUNGNOEN, 2002).

Com as constantes colisões, partículas carregadas penetram num meio até que sua energia

cinética entra em equilíbrio térmico com as partículas do meio, estabelecendo um

“alcance” no meio absorvedor, após um percurso em linha reta ou em zig-zag. As

partículas pesadas, como alfa e fragmentos de fissão, têm uma trajetória praticamente

retilínea dentro do material, ao contrário dos elétrons que têm uma trajetória tortuosa. Para

cada tipo de partícula pode-se definir um alcance, utilizando variações da definição

provenientes de dificuldades experimentais em sua determinação (TAUHATA, 2003).

Alcance, para partículas carregadas, é um conceito mais complexo do que aquele definido

para partículas sem carga, uma vez que estas, em geral, descrevem uma trajetória tortuosa

ao se deslocarem em um meio. Esta complexidade é ainda maior quando se trata de

elétrons de baixa energia, que se deslocam em um pequeno volume até a perda completa de

sua energia. Existe um grande número de estudos que examina a dependência entre energia

e alcance de elétrons com energias desde alguns keV’s até centenas de MeV’s. Porém, a

descrição quantitativa do transporte de elétrons de baixa energia é um problema complexo,

especialmente devido à não disponibilidade de dados confiáveis para as seções de choque

dos diferentes tipos possíveis de interações de elétrons com a matéria

(MEESUNGNOEN, 2002). Os modelos analíticos de transporte de elétrons ou são

incompletos ou apresentam dificuldades até o momento intransponíveis – basta lembrar

que a equação de transporte de Boltzmann é uma equação íntegro-diferencial não linear.

Os códigos de Monte Carlo têm sido a alternativa adotada e são hoje extensivamente

utilizados. Códigos de Monte Carlo de simulação análoga (que simulam o frenamento de

elétrons evento por evento) têm sido uma valiosa ferramenta para explicar, em termos

quantitativos, a ação da radiação no nível molecular (EMFIETZOGLOU, 2003).

O problema tratado neste trabalho consiste na simulação do transporte de elétrons com

energias entre 110 eV e 50 keV, e das partículas secundárias geradas por estes,

provenientes de uma fonte monoenergética isotrópica localizada no interior de um meio

tecido equivalente homogêneo cujas dimensões são muito maiores do que o alcance das

partículas em questão, podendo então ser considerado infinito. E este tem por objetivo

estudar a relação entre a energia inicial e o alcance em meio biológico destes elétrons de

baixa energia. Para tal será utilizado o código de Monte Carlo PENELOPE que simula o

transporte de elétrons e fótons até energias tão baixas quanto 100 eV.

REVISÃO DE LITERATURA

Todo este trabalho fundamenta-se no transporte de radiação – mais precisamente de

elétrons e fótons – através do método de Monte Carlo.

Nas próximas seções será apresentada a teoria relacionada a este assunto. Na seção 1, são

introduzidos os conceitos básicos do método de simulação de Monte Carlo, incluindo uma

breve visão geral dos métodos de amostragem de variáveis aleatórias. Na seção 2 é tratada

a interação da radiação com a matéria. E, finalmente, na seção 3 são apresentadas as

técnicas de Monte Carlo para o caso específico de transporte de radiação (elétrons e

fótons), incluindo a descrição do código PENELOPE.

1 O método Monte Carlo

O método de Monte Carlo é uma técnica para solução numérica de problemas baseada na

simulação de variáveis aleatórias. Sua origem foi marcada pela publicação, em 1949, do

artigo “The Monte Carlo method” e deve sua existência aos esforços de J. von Neumann e

S. Ulam. Através deste método é possível simular tanto processos que dependem de fatores

aleatórios como problemas não ligados a tais fatores, desde que seja possível associar a

estes um ou vários modelos probabilísticos. O campo de aplicação do método é, portanto,

considerável (SOBOL, 1983).

A característica essencial da simulação Monte Carlo é, então, o uso de números e variáveis

aleatórias (SALVAT, 2001). Designa-se variável aleatória aquela que pode assumir

diversos valores, cada um com uma certa probabilidade, não sendo o valor assumido

necessariamente preciso. Mesmo conhecendo-se a probabilidade associada a cada valor

possível de uma variável aleatória não se consegue prever exatamente o valor que esta

assumirá numa experiência, mas se pode prever o caráter do seu comportamento em uma

série de experiências, e a precisão deste será tanto maior quanto mais longa for esta série

(SOBOL, 1983). Esta aleatoriedade pode ser conseqüência da existência de fatores não

controláveis (como ocorre, por exemplo, em jogos de azar) ou em decorrência da natureza

quântica de sistemas e processos microscópicos (por exemplo, desintegração nuclear e

interações da radiação) (SALVAT, 2001).

Apesar do vasto campo de aplicação e da relativa simplicidade deste método, a sua

utilização somente se difundiu após o aparecimento, no final da década de 50, de

computadores capazes de lidar com o considerável volume de cálculos relacionados à

simulação de variáveis aleatórias. Em um computador, estas variáveis são geradas através

de transformações numéricas de números aleatórios (SALVAT, 2001). Nas seções 1.1 e

1.2 serão apresentadas certas propriedades das variáveis aleatórias assim como alguns

métodos de simulação destas.

1.1 Elementos da teoria de probabilidade

Seja x uma variável aleatória contínua que pode assumir qualquer valor no intervalo

maxmin

xxx ≤≤ . A probabilidade de se obter x em um intervalo (a,b) é dado por

{}

NnbxaxP

<< | – razão entre o número n de valores de x que caem neste intervalo

e o número total N de valores gerados de x, no limite

∞

→N . A probabilidade de se obter

x em um intervalo diferencial de comprimento dx em torno de x

pode ser expressa como

, )(} | {

111

dxxpdxxxxxP

( 1.1 )

onde )(xp é a função distribuição de probabilidade (FDP) de x. Uma vez que

probabilidades negativas não têm significado e o valor obtido para x deve estar no intervalo

),(

maxmin

xx , a FDP deve satisfazer as condições:

∫

=≥

max

min

1 )( e 0)(

dxxpxp . ( 1.2 )

Qualquer “função” que obedece estas duas condições pode ser interpretada como uma

FDP. A distribuição uniforme

⎩

⎨

⎧

<<−

≡

maxmin

maxminminmax

ou se 0

se )/(1

)(

maxmin

xxxx

xxxxx

, ( 1.3 )

por exemplo, é freqüentemente utilizada em simulações Monte Carlo. A definição (1.2)

também inclui distribuições peculiares como a de delta de Dirac, )(

xx −

, que é definida

pela propriedade

⎩

⎨

⎧

=−

∫

bx ax

bxa xf

dxxxxf

ouse 0

se )(

)( )(

, ( 1.4 )

para qualquer função )(xf contínua em x

. Uma definição equivalente, e mais intuitiva, é

a seguinte:

),(lim)(

xUxx

xx ∆+∆−

→∆

≡

−

( 1.5 )

que representa a distribuição delta como uma distribuição uniforme centrada no ponto x

no limite em que a largura desta tende a zero. Assim, a distribuição de Dirac descreve uma

variável aleatória que pode assumir um único valor, isto é, uma constante. A FDP de uma

variável aleatória x que assume valores discretos

K x xx ,,

com probabilidades

K p p ,,

pode ser considerada como uma combinação de distribuições delta

∑

−=

xxpxp )( )(

. ( 1.6 )

Distribuições discretas podem então, serem consideradas como formas particulares de

distribuições contínuas.

Dada uma variável aleatória contínua x, a função distribuição cumulativa de x é definida

por:

' )'()(

min

dxxpx

∫

≡

℘

. ( 1.7 )

Esta é uma função não-decrescente de x que varia de

0)(

min

℘

x a 1)(

max

℘

x . No caso

de uma FDP discreta da forma (1.6), )(x

℘

é uma função degrau. Note que a probabilidade

{}

| bxaxP << de obter x no intervalo (a,b) é

∫

℘−℘==<<

abdxxpbxaxP

)()( )( } | { , ( 1.8 )

e que

dxxdxp )()( ℘= .

O n-ésimo momento de )(xp é definido como

∫

max

min

)(

dxxpxx . ( 1.9 )

O momento

x é simplesmente a integral de )(xp , que é igual a um, por definição.

Momentos de ordem maior podem ou não existir O primeiro momento, quando existe, é

chamado de média ou valor esperado da variável aleatória x,

∫

= dxxpxx )( . ( 1.10 )

O valor esperado de uma função )(xf é definido de maneira similar,

∫

≡ dxxpxfxf )( )( )( . ( 1.11 )

Sendo )(xf uma variável aleatória, ela tem sua própria FDP, )( f

, que é tal que a

probabilidade de se ter f em um intervalo de largura df é igual à probabilidade de x estar

no intervalo ou intervalos correspondentes

. Então, se )(xf é uma função

monotonicamente crescente de x (tal que exista uma correspondência um-a-um entre os

valores de x e f ), dffdxxp ) ( )(

= e

)/( )()(

−

= dxdfxpf

. ( 1.12 )

Pode-se mostrar que as definições (1.10) e (1.11) são equivalentes. Se )(xf cresce

monotonicamente com x, a prova é trivial: começando da definição (1.10) e escrevendo

∫

=== dxxpxfdfdfdxxpxfdffff )( )( )/( )( )( ) (

, ( 1.13 )

que concorda com (1.11). Note que o valor esperado é linear, isto é:

)()()()(

22112211

xfaxfaxfaxfa +=+ , ( 1.14 )

onde a

e a

são constantes reais arbitrárias.

Se o primeiro e o segundo momento da FDP )(xp existem, define-se a variância de x (ou

de )(xp ) como

∫

−=−=−≡

222

)()()()var( xxdxxpxxxxx . ( 1.15 )

A raiz quadrada da variância,

[

]

2/1

)var(x≡

, é chamada de desvio padrão, e dá uma

medida da dispersão da variável aleatória, isto é, da largura da FDP. O delta de Dirac é a

única FDP que tem variância nula. De maneira similar, a variância de um função )(xf é

definida como

)()()}(var{ xfxfxf −= . ( 1.16 )

Assim, para uma função constante, axf

)(, af = e

{

}

0var

f .

Variáveis aleatórias bidimensionais

Considerando o caso de uma variável aleatória bidimensional, ),( yx , a FDP (conjunta)

correspondente ),( yxp satisfaz as condições

∫

=≥ 1yxpdy dx yxp ),( e0),( . ( 1.17 )

As FDP’s marginais de x e y são definidas como

Quando f(x) não cresce ou decresce monotonicamente com x, pode haver mais de um valor de x

correspondente a um dado valor de f.

∫

≡≡ dxyxpyq dy yxpxq ),()(e ),()( , ( 1.18 )

isto é, )(xq é a probabilidade de se obter o valor x e qualquer valor y. A FDP conjunta

pode ser escrita como

)|()()|()(),( yxpyqxypxqyxp

, ( 1.19 )

onde

)(

),(

)|(e

)(

),(

)|(

yxp

xyp

yxp

yxp == ( 1.20 )

são as FDP’s condicionais de x e y, respectivamente. Note que )|( yxp é a FDP

normalizada de x para um valor fixo de y.

O valor esperado de uma função ),( yxf é

∫

= ),(),(),( yxpyxfdy dxyxf . ( 1.21 )

Os momentos da FDP são definidos por

∫

= ),( yxpyxdy dxyx

mnmn

. ( 1.22 )

Em particular,

∫

== dxxqxyxpxdy dxx

nnn

)(),( . ( 1.23 )

Novamente, o único momento que é necessariamente definido é

=yx . Quando os

correspondentes momentos existem, as variâncias de x e y são dadas por:

)var(e)var( yyy xxx −=−= . ( 1.24 )

A variância de

),cov(2)var()var()()var(

yxyxyxyxyx ++=+−+=+ , ( 1.25 )

onde

yxxyyx −=),cov( ( 1.26 )

é a covariância de x e y, que pode ser positiva ou negativa. Note que )var(),cov( xxx = e

)var(4)2var( xx = . Quando as variáveis x e y são independentes, isto é, quando

)()(),( ypxpyxp

= , tem-se que

)var()var()var(e0),cov( yxyx yx

. ( 1.27 )

Além disto, para variáveis independentes,

)var()var(}var{

121

yaxayaxa +=+ . ( 1.28 )

(SALVAT, 2001).

1.2 Métodos de simulação de variáveis aleatórias

O primeiro passo para realização de um cálculo de Monte Carlo é a simulação de variáveis

aleatórias com FDP’s específicas. Nesta seção serão descritas diferentes técnicas para a

obtenção de valores aleatórios de uma variável x distribuída no intervalo ),(

maxmin

xx de

acordo com uma dada FDP )(xp ; técnicas que são, inclusive, utilizadas no código

PENELOPE.

1.2.1 Gerador de números aleatórios

Em geral, algoritmos de amostragem randômica são baseados no uso de números aleatórios

uniformemente distribuídos em um intervalo (0,1). Estes números podem ser facilmente

gerados em um computador. Entre os “bons” geradores de números aleatórios disponíveis

atualmente, os mais simples são os chamados geradores congruenciais multiplicativos. Um

exemplo de gerador deste tipo é o seguinte

),12/(),12(mod7

3131

−=−=

− nnnn

R RR

( 1.29 )

que produz uma seqüência de números aleatórios uniformemente distribuídos em (0,1) a

partir de uma dada “semente” )1(

−<

2 R . De fato, a seqüência gerada não é

verdadeiramente aleatória, uma vez que ela é obtida a partir de um algoritmo

determinístico (o termo “pseudo-aleatório” seria mais apropriado), porém é muito pouco

provável que a sutil correlação entre os valores da seqüência tenha um efeito considerável

nos resultados da simulação. O gerador (1.29) é conhecido por ter boas propriedades

randômicas, porém, a seqüência gerada por este é periódica, com um período da ordem de

que, com os recursos computacionais atuais, não é grande o suficiente para evitar a

reinicialização dos números aleatórios em uma simulação. Se o período não for

suficientemente grande para evitar esta reinicialização, é necessário o uso de algoritmos

mais sofisticados do que os congruenciais. O PENELOPE utiliza o gerador de L’Ecuyer

(1988) (citado por SALVAT, 2001)implementado no FORTRAN77, que possui um

período da ordem de 10

, virtualmente infinito na prática

1.2.2 Método da transformada inversa

A função distribuição cumulativa de )(xp , equação. (1.7), é uma função não-decrescente

de x e, portanto, possui uma função inversa

)(

−

℘

. A transformação )(x℘=

define uma

nova variável aleatória que assume valores no intervalo (0,1) (FIGURA 1.1) Devido à

correspondência entre os valores de x e

, as FDP’s de

, )(

p , e x, )(xp , são

relacionadas por dxxpdp )( ) ( =

. Assim,

)(

)()()(

⎟

⎠

⎞

⎜

⎝

⎛

℘

⎟

⎠

⎞

⎜

⎝

⎛

−−

xpp

, ( 1.30 )

ou seja,

é uniformemente distribuído no intervalo (0,1).

Então, sendo

um número aleatório, a variável x definida por

)(

−

℘=

x é aleatoriamente

distribuída no intervalo ),(

maxmin

xx com FDP )(xp (FIGURA 1.1). Isto estabelece um

método prático de gerar valores aleatórios de x utilizando um gerador de números

aleatórios uniformemente distribuído em (0,1). A aleatoriedade de x é garantida por aquela

. Note que x é a raiz (única) da equação

∫

dxxp

min

' )'(

, ( 1.31 )

que será chamada de equação de amostragem da variável x. Este procedimento de

amostragem aleatória é conhecido como método da transformada inversa, e é adequado

para FDP’s dadas por expressões analíticas simples tais que a equação de amostragem

(1.31) possa ser resolvida analiticamente.

FIGURA 1.1: Amostragem aleatória a partir de uma distribuição p(x)

usando o método da transformada inversa.

Transformada inversa numérica

O método da transformada inversa pode ser eficientemente utilizado na amostragem

aleatória a partir de distribuições contínuas )(xp dadas na forma numérica ou cuja

equação de amostragem correspondente seja muito complicada para ser resolvida

analiticamente. Para aplicar este método, os valores da função distribuição cumulativa

)(x

℘

devem ser calculados em alguns pontos x

do domínio de )(xp e armazenados na

memória do computador. A equação de amostragem

℘

)(

x pode então ser resolvida

por interpolação na tabela ),(

, onde )(

℘

≡

, tratando-se

como variável

independente.

Distribuições discretas

O método da transformada inversa pode também ser aplicado a distribuições discretas.

Considerando uma variável aleatória que pode assumir os valores discretos N , x

K,1

com probabilidades

p p ,,

K , respectivamente, a FDP correspondente pode ser escrita

como

∑

−=

ixpxp

)( )(

, ( 1.32 )

onde )(x

é a distribuição de Dirac. Aqui, assume-se que )(xp é definida em um

intervalo ),( ba com

1<a e Nb > . A função distribuição cumulativa correspondente é

[

]

∑

℘

)( , ( 1.33 )

onde

[]

x denota a parte inteira de x. Note que 0)(

℘

x quando 1<x . Assim, equação

(1.31) leva à fórmula de amostragem:

pp j x

ppp x

p x

∑∑

−

=≤<=

+≤<=

≤

211

1se

se2

se1

( 1.34 )

Então, para obter x, basta gerar um número aleatório

e escolher x igual ao índice i tal

que

≤

, ( 1.35 )

onde

∑

==+===

pP ..., ppP pP P

1213121

1,,,0 . ( 1.36 )

(SALVAT, 2001).

1.2.3 Método de Walker

Em 1977, Walker definiu o mais eficiente dos métodos de amostragem a partir de

distribuições discretas, que sorteia um valor com apenas uma comparação. A idéia básica

do método de Walker pode ser facilmente compreendida recorrendo-se a argumentos

gráficos (SALVAT, 1987 citado por SALVAT, 2001). Para tal, a FDP (1.32) será

representada como um histograma construído com

N colunas de largura N1 e altura

(FIGURA 1.2). As barras do histograma podem ser cortadas convenientemente de maneira

que as peças resultantes sejam encaixadas formando um quadrado de lado unitário tal que

cada linha vertical cruze, pelo menos, duas peças diferentes. Este arranjo pode ser

construído sistematicamente selecionando-se a menor e a maior barra do histograma (aqui,

identificadas com os índices l e j, respectivamente) e cortando-se a maior para completar a

menor, que é mantida inalterada. Para registrar a transformação realizada, rotula-se a peça

adicionada com o valor jK

= , que identifica sua posição original no histograma, e

FIGURA 1.2: Representação gráfica do método da transformada inversa (topo) e do método de

Walker para amostragem aleatória a partir de uma distribuição discreta. Neste exemplo, a variável

aleatória pode assumir os valores i = 1, 2, 3 e 4 com probabilidades relativas 1, 2, 5 e 8.

introduz-se o valor “de corte”

F definido como a altura da menor peça na l-ésima coluna

do quadrado resultante; esta peça mantém o rótulo l. Evidentemente, a repetição deste

processo completará o quadrado (após

−

N passos). Note que as probabilidades

podem ser reconstruídas a partir dos rótulos e valores de corte:

),( )1(

jii

KiFFNp

∑

≠

−+=

, ( 1.37 )

onde ),( ji

é o delta de Kronecker ( 1

e 0

caso contrário). A amostragem de

uma variável aleatória

através do método de Walker ocorre como a seguir: Gera-se dois

números aleatórios,

, definindo um ponto aleatório ),(

, que é uniformemente

distribuído no quadrado. Se

),(

estiver em uma peça rotulada com o índice i, define-se

= . Obviamente, a probabilidade de se obter i como resultado é igual à fração da área do

quadrado ocupada pelas peças de rótulo i, que coincide com

p .

Como descrito acima, o algoritmo de Walker requer a geração de dois números aleatórios

para cada valor sorteado para

. Mas é possível gerar os valores de

a partir de um único

número aleatório: assumindo que as

N colunas do quadrado são alinhadas

consecutivamente formando um segmento de comprimento

N (parte inferior da

FIGURA 1.2). Desta forma, para sortear um valor para

, basta gerar um número

aleatório N

, que é uniformemente distribuído em ),0( N e determinar uma das peças do

segmento. O resultado da amostragem é igual ao rótulo da peça selecionada.

Resumidamente, o algoritmo procede como a seguir:

(i) Gera-se um número aleatório

e define-se 1 +

(ii) Faz-se

[]

Ri = e iRr

−

(iii) Se

Fr > ,

Kx = .

(iv) Caso contrário, i

Desta forma, a amostragem de

envolve a geração de apenas um número aleatório e uma

comparação (independente do número

N de resultados possíveis). O preço que se paga

por esta simplificação é o aumento da quantidade de memória necessária para executar o

algoritmo: ao invés de um único banco de dado,

p (ou

P ), são usados dois,

K e

F .

Infelizmente o cálculo de rótulos e valores de corte é bastante complexo, o que limita a

aplicabilidade do algoritmo de Walker para distribuições que se mantêm constantes

durante a simulação (SALVAT, 2001).

1.2.4 Métodos de rejeição

Os algoritmos de rejeição baseiam-se na simulação de uma variável aleatória x – a partir de

uma certa distribuição diferente de )(xp – e em um teste que determina se esta variável

será rejeitada ou não. A compreensão destes métodos pode ser feita através de simples

argumentos gráficos (FIGURA 1.3): Utilizando-se o método da transformada inversa ou

qualquer outro método de amostragem aplicável, gera-se valores aleatórios para a variável

x a partir de uma FDP )(x

. Para cada valor de x, gera-se um valor aleatório para y

uniformemente distribuído no intervalo ))(,0( xC

, onde C é uma constante positiva.

Desta forma, os pontos ),( yx gerados estarão uniformemente distribuídos na região A do

plano delimitada pelo eixo x (0

y ) e pela curva )(xCy

. De maneira recíproca, se

forem gerados, de alguma forma, pontos ),( yx aleatórios uniformemente distribuídos em

A, a coordenada x destes pontos será uma variável aleatória distribuída de acordo com

)(x

(independente do valor de C). Considerando que a distribuição )(x

é tal que

)()( xpxC ≥

para algum valor de 0>C e que os pontos aleatórios ),( yx estão

distribuídos uniformemente em A como descrito acima, se os pontos em que )(xpy >

forem rejeitados, os pontos aceitos (com )(xpy

≤

) estarão uniformemente distribuídos na

região entre o eixo x e a curva )(xpy

e, então, as coordenadas x correspondentes estarão

distribuídas de acordo com )(xp .

Um método de rejeição pode então ser representado pela FDP )(xp escrita como

)( )( )( xrxCxp

, ( 1.38 )

FIGURA 1.3: Amostragem aleatória a partir de uma distribuição p(x) utilizando um método de rejeição.

onde )(x

é uma FDP escolhida tal que facilite a simulação de x (utilizando, por exemplo,

o método da transformada inversa), C é uma constante positiva e a função )(xr satisfaz a

condição 1)(0 ≤< xr .

Com os argumentos gráficos utilizados acima, fica claro que este algoritmo fornece valores

de x distribuídos de acordo com )(xp . A eficiência do algoritmo, isto é, a probabilidade de

se aceitar um valor gerado para x, é:

∫

dxxxr

)( )(

πε

. ( 1.39 )

Graficamente, a eficiência é igual à razão entre as áreas sob as curvas )(xpy

)(xCy

= , que são iguais a 1 e C, respectivamente. Como 1)(

≤

xr , para uma dada )(x

a constante C deve satisfazer a condição )()( xpxC ≥

para qualquer x. O menor valor

possível de C, que ocorre quando )()( xpxC

, fornece a melhor eficiência.

A FDP )(x

deve ser escolhida de maneira tal que o algoritmo de amostragem seja o mais

rápido possível. Uma grande eficiência também é desejável, mas não é decisiva. Cem por

cento de eficiência é obtida somente quando )()( xpx

, qualquer outra FDP resultará em

menor eficiência (porém a amostragem a partir desta FDP é justamente o problema que se

deseja resolver). A utilidade do método de rejeição está no fato de que a perda na

eficiência pode ser compensada com a simplicidade da simulação de x a partir de )(x

invés de )(xp . A desvantagem deste método é a necessidade de vários números aleatórios

ξ para gerar cada valor de x (SALVAT, 2001).

1.2.5 Métodos de composição. Variáveis bidimensionais

Como visto anteriormente na seção 1.1 [equações. (1.18), (1.19) e (1.20)], para uma

variável aleatória bidimensional, a função distribuição de probabilidade pode ser escrita

como

)|()(),( yxpyqyxp

, ( 1.40 )

onde

∫

=≡

)(

),(

)|(e),()(

yxp

yxp dxyxpyq ( 1.41 )

são, respectivamente, a FDP marginal de y e a FDP condicional de x.

Para gerar pontos ),( yx aleatórios a partir de ),( yxp basta obter y de )(yq e, então, x de

)|( yxp . Portanto, variáveis aleatórias bidimensionais podem ser geradas usando-se

métodos de amostragem de variáveis simples. Isto vale também para variáveis n-

dimensionais, uma vez que uma FDP n-dimensional pode sempre ser expressa como o

produto entre uma distribuição marginal de uma variável simples e uma FDP de dimensão

)1( −n .

Da definição da FDP marginal de x,

∫

=≡ dyyxpyqdyyxpxq )|( )( ),()( , ( 1.42 )

fica claro que, determinando y a partir de )(yq e, então, x de )|( yxp , os valores gerados

para x são distribuídos de acordo com )(xq . Esta é a base dos métodos de composição, que

são aplicáveis quando a distribuição, )(xp , a ser simulada é uma mistura de várias FDP’s.

Esta técnica é muito útil na simulação de valores aleatórios a partir de distribuições

complexas, que podem ser expressas como uma combinação de distribuições mais simples

que são, por sua vez, facilmente tratadas pelo método da transformada inversa ou por

algum método de rejeição.

Para distribuições analiticamente simples, cuja função distribuição cumulativa possui

inversa, o método da transformada inversa é, em geral, satisfatório. Este é o caso de

algumas poucas distribuições elementares (por exemplo, as distribuições uniforme e

exponencial). Este método é também adequado para distribuições discretas e para FDP’s

contínuas dadas na forma numérica. Combinando-se os métodos da transformada inversa,

de rejeição e de composição, é possível desenvolver algoritmos de amostragem para

praticamente todas FDP’s (SALVAT, 2001).

1.3 Integração Monte Carlo

Todo cálculo de Monte Carlo, pelo menos no sentido formal, é equivalente à integrações

(JAMES, 1980 citado por SALVAT, 2001). Esta equivalência permite estabelecer

formalmente uma teoria para as técnicas de Monte Carlo. Um importante aspecto em

simulações é a determinação das incertezas estatísticas das grandezas calculadas. Nesta

seção serão estabelecidas as fórmulas básicas para tais cálculos considerando-se o cálculo

de Monte Carlo mais simples: o de uma integral unidimensional. Evidentemente estes

resultados também serão válidos para integrais de dimensões maiores.

A integral

∫

dxxFI

)( , ( 1.43 )

pode ser escrita na forma de um valor esperado

∫

≡= fdxxpxfI )( )( , ( 1.44 )

introduzindo uma FDP arbitrária )(xp e fazendo

)()()( xpxFxf

(assume-se que

0)( >xp em ),( ba e 0)( =xp fora deste intervalo). O cálculo Monte Carlo da integral I é

muito simples: gera-se um grande número N de pontos aleatórios x

a partir da FDP )(xp

acumulando, simultaneamente, a soma dos valores )(

xf em um contador. No final do

cálculo o valor esperado de f será estimado como

∑

≡

)(

. ( 1.45 )

A lei dos grandes números diz que, quando N é suficientemente grande,

ade)probabilid (em If → . ( 1.46 )

Em termos estatísticos, isto significa que

f , o resultado de Monte Carlo, é um avaliador

consistente para a integral (1.43). Isto é válido para qualquer função )(xf finita e com um

número finito de descontinuidades.

A lei dos grandes números (1.46) pode ser restabelecida como

∑

∞→

)(

lim . ( 1.47 )

Aplicando esta lei à integral que define a variância de )(xf [equação. (1.16)]

∫

−=

)( )( } )( var{ fdxxpxfxf , ( 1.48 )

tem-se que

[]

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−=

∑∑

∞→

)(

lim)}(var{ . ( 1.49 )

A expressão entre as chaves é uma estimativa para variância de )(xf .

É visível que cálculos de Monte Carlo distintos (com seqüências de números aleatórios x

diferentes e independentes) resultarão em diferentes valores para

f . Isto implica que os

resultados de um código Monte Carlo é afetado por incertezas estatísticas, similares

àquelas encontradas em experimentos laboratoriais, que precisam ser calculadas

adequadamente para determinar a “precisão” do resultado de Monte Carlo. Para tal,

considera-se

f como uma variável aleatória cuja FDP é, em princípio, desconhecida. Isto

significa que a média e a variância de

f são dadas por

∑∑

===

)(

( 1.50 )

{} {}

∑∑

⎥

⎦

⎤

⎢

⎣

⎡

)(var

)(

var)var( , ( 1.51 )

onde foram utilizadas as propriedades (1.14) e (1.28) do valor esperado e da variância,

respectivamente. O desvio padrão de

f ,

{

}

)(var

)var( =≡

, ( 1.52 )

fornece o valor da incerteza estatística da estimativa Monte Carlo de

f . O resultado (1.52)

tem uma importante implicação prática: para reduzir a incerteza estatística por um fator de

10, é necessário aumentar o tamanho da amostra N por um fator de 100, o que determina

um limite na precisão que pode ser alcançada em um cálculo em função da limitação

computacional.

Recorrendo ao teorema central do limite (SOBOL, 1983), tem-se que, no limite

∞

→N , a

FDP

)( fp de f , é uma distribuição normal (Gaussiana) com média f e desvio padrão

⎟

⎠

⎞

⎜

⎝

⎛

−

−=

)(

exp

)(

πσ

. ( 1.53 )

Das propriedades da distribuição normal, segue que, para valores de N grandes o

suficiente, para os quais o teorema se aplica, o intervalo

± contém o valor exato

f – isto é, o valor da integral I – com probabilidade de 68,3% se 1=n , de 95,4% se

2=n e de 99,7% se 3=n . Na regra utiliza-se 3

n , uma vez que a diferença entre a

probabilidade 99,7% e 100% é pouco significativa e, portanto, é quase certo que, em uma

experiência, o desvio entre

f e f seja inferior a

3 .

O teorema central do limite é uma poderosa ferramenta uma vez que ele prevê que os

valores de f seguem uma distribuição específica; porém ele se aplica apenas

assintoticamente. O número mínimo N de amostras necessárias para se aplicar o teorema

seguramente depende do problema em questão. Se o terceiro momento central de f,

∫

−≡ dxxpfxf )(])([

, ( 1.54 )

existe, o teorema é satisfeito quando

σµ

<< . ( 1.55 )

Em geral, aconselha-se a análise da distribuição de

f para certificar-se da aplicabilidade

do teorema central do limite. Porém, na maioria dos cálculos de Monte Carlo, erros

estatísticos são estimados simplesmente assumindo que o teorema é satisfeito,

independente do tamanho da amostra.

Cada )(xp possível define um algoritmo de Monte Carlo para calcular a integral I,

equação (1.43). O mais simples deles é obtido utilizando-se a distribuição uniforme

)(1)( abxp −= . Evidentemente, )(xp determina não somente a densidade dos pontos

da amostra, mas também a magnitude da variância

{

}

)(var xf , equação (1.48),

{}

∫∫

⎥

⎦

⎤

⎢

⎣

⎡

−=−

⎥

⎦

⎤

⎢

⎣

⎡

dxI

xFIdx

xpxf

)(

)( )( var . ( 1.56 )

Como medida da eficácia de um algoritmo de Monte Carlo é comum o uso da eficiência

definida por

[]

= , ( 1.57 )

onde T é o tempo computacional para realização da simulação (ou qualquer outra medida

do esforço computacional). Uma vez que

e T são aproximadamente proporcionais a

1−

N e N respectivamente,

é uma constante (independente de N) na média.

Os chamados métodos de redução de variância são técnicas que ajudam a otimizar a

eficiência da simulação através da escolha adequada da FDP )(xp . Porém, uma redução na

variância não necessariamente implica em uma maior eficiência; se um algoritmo de

Monte Carlo baseado em uma certa FDP )(xp tem variância menor do que aquela de uma

distribuição uniforme, mas leva mais tempo para gerar os valores de x a partir de )(xp ,

este algoritmo de “variância reduzida” pode ser menos eficiente do que aquela que

considera a distribuição uniforme. Portanto, deve-se evitar o uso de FDP’s difíceis de lidar

na amostragem O aumento da eficiência de algoritmos é uma parte importante e delicada

da técnica de Monte Carlo. (SALVAT, 2001).

2 Interação da radiação com a matéria

Nesta seção são apresentados os principais mecanismos de interação de fótons e elétrons

com a matéria relativos à faixa de energia considerada neste trabalho. No apêndice B são

descritos alguns dos termos, grandezas e unidades relacionados ao assunto para que o leitor

possa recorrer se necessário. No apêndice C é feita uma revisão sobre cinemática de

colisão.

2.1 Interação de fótons com a matéria

As interações primárias de fótons com a matéria envolvem a produção de partículas

secundárias carregadas, geralmente elétrons. É predominantemente através da interação

destas partículas que ocorre a transferência da energia para o material absorvedor,

usualmente através da ionização dos átomos deste; por isso, os efeitos biológicos

provenientes da interação de fótons são na verdade os efeitos biológicos provocados pela

interação de partículas secundárias com o meio.

Um fóton pode interagir com a matéria por qualquer um dos vários mecanismos

competitivos alternativos. A interação pode ser com o átomo como um todo (efeito

fotoelétrico e espalhamento Rayleigh), com um elétron do átomo (efeito Compton e

produção de pares no campo eletromagnético de um elétron) ou com o núcleo atômico

(produção de pares, espalhamento elástico, foto-desintegração e produção de mésons). A

probabilidade de ocorrência de cada um destes processos pode ser expressa como uma

seção de choque de colisão por átomo, por elétron ou por núcleo no absorvedor. A soma de

todas estas seções de choque, normalizada por átomo, é então a probabilidade de um fóton

incidente interagir de alguma forma enquanto passa através de um absorvedor muito

estreito que contém um átomo por centímetro quadrado de área normal à direção de

incidência do fóton.

A seção de choque de colisão total por átomo multiplicada pelo número de átomos por

centímetro cúbico do absorvedor é igual ao coeficiente de atenuação linear

por

centímetro percorrido no absorvedor. A fração de fótons incidentes que passa através de

um absorvedor de largura

e densidade

sem sofrer qualquer tipo de interação é

[]

−exp ou

[]

))(/(exp x

− , onde

/ é o coeficiente de atenuação por massa. A

atenuação pode ocorrer por algum processo puramente elástico, como o espalhamento

Rayleigh ou espalhamento por ressonância nuclear (Mossbauer), onde o fóton é

simplesmente defletido e não perde energia para o meio. Neste caso, estaria envolvido no

processo somente um coeficiente de espalhamento,

. Por outro lado, em uma interação

fotoelétrica, toda energia do fóton incidente é absorvida por um átomo do meio, não

havendo fóton espalhado residual; aqui, a atenuação da radiação primária é devida à

absorção completa da energia do fóton incidente. O caso intermediário, e de grande

importância, é o efeito Compton em que parte da energia é absorvida e aparece no meio

como energia cinética de recuo de um elétron, enquanto o restante da energia incidente não

absorvida se apresenta como um fóton Compton espalhado. Portanto, o espalhamento

envolve a deflexão do fóton incidente e a absorção envolve a conversão da energia do

fóton incidente em energia cinética de partículas carregadas, que posteriormente é

dissipada através de colisões ao longo de suas trajetórias. O coeficiente de absorção

tem um conceito mais restrito do que o coeficiente de atenuação

, que é a soma do

coeficiente de espalhamento

e do próprio coeficiente de absorção

. Em muitos casos

práticos,

é simplesmente o coeficiente de espalhamento Compton

(ATTIX, 1968).

No domínio das energias encontradas mais freqüentemente (0,01 a 10 MeV) a maioria dos

efeitos de interação são devidos a somente três dos vários processos competitivos: o efeito

Compton, o efeito fotoelétrico e a produção de pares (ATTIX, 1968). As características

essenciais destes processos são (ROSSI, 1994):

a. Efeito fotoelétrico (EF) e produção de pares (PP) são eventos de absorção, ou

seja, o fóton desaparece após estas interações. No espalhamento Compton (C),

como o próprio nome define, o fóton sofre apenas um espalhamento elástico.

b. A dependência das seções de choque,

, com o número atômico do material

alvo é:

254

~ ,~ ,~ ZZZZ

PPCEF

σσσ

− .

c. Quando se trata de tecido biológico (

Z ), as interações fotoelétrica,

Compton e produção de pares são dominantes, respectivamente, nos intervalos

de energia: < 50 keV, (50 keV,10 MeV) e > 10 MeV (FIGURA 2.1).

A FIGURA 2.1 apresenta a importância relativa destas três principais formas de interação

em função da energia

h do fóton incidente e do número atômico

do material

atenuador. Para qualquer

o efeito Compton predomina para energias do fóton entre 0,8 e

4 MeV; para materiais de baixo

o efeito Compton predomina sobre um domínio maior

de energia. Para Z moderadamente grande, a interação fotoelétrica é dominante para

h pequeno e a produção de pares para

h grande (ATTIX, 1968)

No caso deste trabalho os mecanismos de interação relevantes são o efeito fotoelétrico e

espalhamento Compton que são descritos a seguir:

2.1.1 Efeito fotoelétrico

Neste processo um fóton de energia

E colide com um elétron que orbita um átomo do

material. O fóton é absorvido e o elétron atômico é liberado com energia cinética

WET

−

, ( 2.1 )

onde

W é a energia de ligação da camada da qual o elétron foi removido. Para estimar

aproximadamente o valor de

W pode-se usar as seguintes fórmulas empíricas (em

unidades de 13.5 eV) (FERMI, 1950 citado por ROSSI, 1994):

9 M camada

4 L camada

K camada

13)(Z

5)(Z

1)(Z

−

( 2.2 )

A seção de choque fotoelétrica,

, é aproximadamente proporcional a

EZ ; porém, à

medida que a energia do fóton aumenta ela aproxima do valor em que o balanço de energia

[equação. (2.1)] permite a remoção de um elétron de uma camada mais interna, ocorrendo

um aumento característico na seção de choque chamado de limite K, L, M, ... . Este

aumento repentino em

é devido ao acréscimo no número de elétrons disponíveis e a

FIGURA 2.1: Importância relativa dos três principais tipos de interação de fótons com a

matéria em função da energia h

do fóton.e do número atômico Z do material.

uma propriedade comum deste processo em que quanto mais fortemente ligados os

elétrons, maior a seção de choque.

A distribuição angular dos elétrons ejetados (fotoelétrons) depende de

E . Em baixas

energias a maioria dos fotoelétrons é emitida ao longo da direção do vetor campo elétrico

do fóton, isto é, a 90º da direção de incidência do feixe. À medida que

E aumenta o

fotoelétron é progressivamente emitido em ângulos menores (tendendo à direção do feixe

incidente).

A vacância deixada pela emissão de um fotoelétron é preenchida por um elétron de uma

camada mais externa e a energia originada nesta transição é liberada através da emissão de

raios-X (fluorescência) ou pela produção de um elétron Auger das camadas L ou M. A

razão entre as probabilidades destes dois processos é chamada “fluorescence yield” e

cresce monotonicamente com Z (ROSSI, 1994).

2.1.2 Espalhamento Compton

Espalhamento Compton é uma colisão elástica entre um fóton e um elétron livre em

repouso. Para elétrons atômicos usualmente assume-se que a energia do fóton é muito

maior do que a energia de ligação, w, dos elétrons atômicos e, por questões práticas, eles

podem ser considerados livres e em repouso.

A geometria do espalhamento Compton é mostrada na FIGURA 2.2. A energia do fóton

espalhado, '

E , é dada por:

)cos1(1

θα

−+

( 2.3 )

onde

= . ( 2.4 )

FIGURA 2.2: No efeito Compton um fóton (de momento cE

) é espalhado por um elétron

atômico. O elétron e o fóton são espalhados com momentos

p e cE '

, respectivamente.

O ângulo de espalhamento do elétron é dado por:

tgtg

= , ( 2.5 )

e, portanto, o elétron é sempre espalhado na direção do fóton incidente. A energia do

elétron,

T , é dada por:

γγ

EET

−

( 2.6 )

)cos1(1

)cos1(

θα

−+

−

= ET

. ( 2.7 )

A energia do elétron é máxima quando o ângulo de espalhamento do fóton é

º180=

)(

max

= ET

. ( 2.8 )

A seção de choque diferencial

)(cos

dd , para o espalhamento de um fóton na direção

é dada pela fórmula de Klein-Nishima (ROSSI, 1952 citado por ROSSI, 1994):

⎭

⎬

⎫

⎩

⎨

⎧

−++

−

−+

)]cos1(1)[cos1(

)cos1(

)]cos1(1[

cos1

)(cos

θαθ

θα

( 2.9 )

onde

mcer

= é o raio clássico do elétron.

A equação abaixo descreve a seção de choque Compton total:

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−

)21(

)21ln()21ln(

)1(21

πσ

r . ( 2.10 )

(ROSSI, 1994).

2.2 Interação de partículas carregadas com a matéria

As partículas carregadas podem ser classificadas como elétrons e partículas carregadas

pesadas. A diferença fundamental entre elétrons e todas as outras partículas carregadas está

na massa de repouso, que no caso do elétron (

m ) é muito pequena.

Ao atravessar um meio, as partículas carregadas perdem energia principalmente através de

interações Colombianas com os elétrons do material. Elétrons incidentes colidindo com

elétrons podem perder uma grande fração de sua energia em uma única interação e

portanto podem ter grandes ângulos de deflexão; partículas pesadas incidentes perdem

apenas pequena fração de sua energia, tipicamente no máximo

Mm4 , onde

é a massa

da partícula pesada, e serão defletidas de ângulos muito pequenos.

Para ambos os tipos de partículas as principais interações com a matéria se dividem

basicamente em três grupos: interações com os elétrons individuais de átomos e moléculas,

interação com o núcleo e interações com o átomo como um todo. A última ocorre para

elétrons com energias abaixo das energias de excitação dos átomos e para partículas

pesadas com velocidades muito pequenas,

(

)

zcv 03.0

≤

, onde z é o número atômico

do íon incidente. Neste caso, a partícula interage com o sistema acoplado do núcleo com os

elétrons atômicos

, em contraste com a situação de altas velocidades em que o

acoplamento entre os elétrons e o resto do átomo tem efeito desprezível na transferência de

energia. As interações com o núcleo incluem reações nucleares, espalhamento nuclear e

espalhamento de Rutherford, todos estes contribuindo para o fenômeno de espalhamento

múltiplo

Estes processos consistem de interações individuais que possuem uma distribuição para as

possíveis transferências de energia. Por esta razão, após a passagem de um grupo de

partículas monoenergéticas através de um material de uma dada espessura, elas irão

apresentar uma dispersão de suas energias (“energy straggling”). Se a espessura é grande o

suficiente, algumas das partículas irão atravessar o material e outras não (“range

straggling”).

Para energias consideravelmente maiores do que a energia associada à massa de repouso

da partícula, a produção de radiação de bremsstrahlung e Cerenkov

tem grande relevância.

A geração de fótons de bremsstrahlung ocorre quando partículas carregadas são

desaceleradas, principalmente no campo Colombiano de um núcleo. Este processo pode

ocorrer também com partículas de baixa energia, porém a contribuição deste processo na

perda de energia de partículas de baixa energia é, em geral, tão pequena que pode ser

desprezada. Para elétrons, a perda de energia por bremsstrahlung começa a ser relevante

para energias maiores que 10 MeV em chumbo e por volta de 100 MeV em água.

Em geral, estas colisões são inadequadamente chamadas de colisões nucleares.

O efeito de várias pequenas deflexões é chamado espalhamento múltiplo. Elétrons sofrem muito mais

espalhamento múltiplo do que partículas pesadas.

A radiação Cerenkov é uma espécie de “onda de choque” eletromagnética, na região do visível, que ocorre

quando uma partícula carregada atravessa uma substância com velocidade maior do que a velocidade da luz

neste material.

Durante o frenamento de partículas carregadas em qualquer material, alguns átomos (ou

moléculas) do meio são excitados, isto é, elétrons são transferidos do estado fundamental

para estados excitados. Alguns dos átomos podem também ser ionizados; isto significa que

o elétron é completamente separado do átomo, deixando para trás um íon positivo.

Geralmente, este elétron secundário irá adquirir uma quantidade

de energia cinética. A

energia de recuo do íon é relativamente pequena e é desprezada na maioria dos casos. Se

for suficientemente grande para permitir a observação do elétron secundário através de

suas próprias interações, este elétron é chamado de raio

. Muitas vezes, o elétron

secundário se liga a uma molécula neutra, formando um íon negativo. Isto ocorre com

grande probabilidade em alguns gases como oxigênio ou vapor d’água, e com pequena

probabilidade em outros, como argônio ou nitrogênio. Um íon negativo (ou elétron livre) e

o íon positivo correspondente são chamados de par de íons. A energia média necessária

para a produção de um par de íons é normalmente chamada

W (ATTIX, 1968).

2.2.1 Poder de frenamento

Em função dos vários processos que contribuem para o frenamento de partículas

carregadas, é conveniente tratar o poder de frenamento em três intervalos de energias

distintos:

a. Baixas energias, (i) abaixo de

21 z MeV para partículas capazes de “atingir”

elétrons orbitais, e (ii) abaixo de aproximadamente

10 Mc

−

para os demais.

b. Energias intermediárias.

c. Altas energias, energia cinética bem acima de

Mc para todas as partículas.

Esta divisão em intervalos de energia permite um tratamento uniforme ao desprezar ou

agrupar certos efeitos.

Em baixas energias, os íons positivos irão capturar e perder elétrons, ficando

completamente ionizados apenas parte do tempo. Uma partícula de carga

z é capaz de

produzir uma ionização completa (arrancando todos elétrons) somente quando possui uma

energia de aproximadamente

21 z MeV. Colisões nucleares (colisões elásticas entre

partícula e o átomo como um todo) são importantes para todas partículas cuja energia é

menor que

10 Mc

−

e requerem um tratamento diferente daquele das colisões inelásticas

com elétron.

Na região de energias intermediárias, a principal interação é a que ocorre entre a partícula

incidente e os elétrons do meio. Os elétrons são elevados para estados de maior energia

(excitação) ou para estados não ligados. O processo de colisão é bem compreendido e, a

princípio, o poder de frenamento (perda de energia média por unidade de comprimento da

trajetória) pode ser calculado teoricamente. As interações nucleares ocorrem com uma

freqüência relativamente pequena (da ordem de 1 a cada 10

interações com elétrons), mas

no caso das partículas carregadas pesadas seus efeitos são relevantes para as energias

maiores da região intermediária. Reações nucleares começam a ser mais e mais

importantes à medida que a energia aumenta. Desta forma, o frenamento é melhor descrito

como um certo número de eventos catastróficos, muito mais parecido com o processo de

absorção de fótons ou nêutrons do que com o processo de frenamento relativamente

contínuo. Praticamente não há sentido em se falar de alcance para um próton de

1000 MeV.

Para energias muito maiores do que

Mc , a emissão de bremsstrahlung é o processo de

maior relevância. Além disto, nestas energias, as principais propriedades do material

(constante dielétrica, etc.) começam a influenciar fortemente o processo de colisão

proporcionando a redução do poder de frenamento. Este fenômeno é chamado de efeito de

densidade ou efeito de polarização.

Fórmula fundamental para o poder de frenamento

Em energias intermediárias, o principal mecanismo de interação de partículas carregadas

com o meio é a excitação e ionização de elétrons do material. Este é um problema

fundamentalmente de mecânica quântica, porém, os resultados obtidos através de um

tratamento clássico que usa órbitas de Bohr para os elétrons são válidos em uma grande

variedade de experimentos de penetração de partículas na matéria. Este tratamento clássico

tem sido usado para todos tipos de partículas carregadas, mas no caso de elétrons são

necessárias algumas modificações no resultado final.

FIGURA 2.3: Colisão de uma partícula positiva pesada (massa M) com um elétron (massa m).

Considere o problema: Uma partícula incidente de carga

ze, massa

, energia cinética

, e velocidade v , movendo-se ao longo do eixo

de um sistema de coordenadas

cartesianas, colide com um elétron livre de carga

−

e massa m situado no eixo y , como

na FIGURA 2.3. Se

for muito maior que m , a deflexão da partícula incidente pode ser

desprezada. Inicialmente, será assumido que o elétron não se movimenta durante a colisão.

O parâmetro de impacto,

b , é definido como sendo a menor distância entre as duas

partículas na ausência de uma interação. De acordo com a segunda lei de Newton, a força

Colombiana

(

)

(

)

22222

cos bzerzeF

== irá transferir um momento q para o elétron:

∫

= dt Fq ( 2.11 )

Com as suposições feitas, a componente de

F na direção x não contribui para a integral, e

pode-se escrever:

∫

== dt

mvq

cos

( 2.12 )

uma vez que

cos F=

F . Para substituir dt na integral, é usada a derivada de

)( bxtg =

em relação ao tempo, isto é, dtbxddttgd /)(/)( =

, onde

)/)(cos/1(/)(

dtddttgd

θθθ

= e bvdtdxbdtbxd /)/)(/1(/)(

. Substituindo as duas

últimas equações na primeira:

bvdtd /)/)(cos/1(

θθ

dvbdt cos)/(

2−

= .

Substituindo-se então

dt na equação (2.12):

º90

cos ==

∫

−

θθ

. ( 2.13 )

Daí a energia

T transferida para o elétron, nesta aproximação clássica, é:

vmb

== . ( 2.14 )

Note que

T depende apenas da velocidade v e da carga z da partícula incidente e não de sua

massa

M. Usando )/(

mcer = (o “raio clássico do elétron”) e introduzindo a energia

cinética

E da partícula incidente tem-se que

⎟

⎠

⎞

⎜

⎝

⎛

. ( 2.15 )

A energia cinética

T transferida para um elétron, em termos da energia de repouso deste, é

então inversamente proporcional à energia cinética

E da partícula incidente (em termos da

sua massa de repouso), inversamente proporcional ao quadrado do parâmetro de impacto

e proporcional ao quadrado da carga

z da partícula.

Existe um limite para a quantidade de energia que pode ser transferida para o elétron. A

velocidade máxima que um elétron pode adquirir é

vu 2

max

, ( 2.16 )

ou seja, duas vezes a velocidade da partícula incidente. A energia máxima transferida é

então

maxmax

mvmuT == . ( 2.17 )

A equação (2.14) na verdade não é válida para os pequenos valores de

b que resultariam

nesta transferência de energia

max

, mas a conclusão a respeito de T

max

ainda é correta. Da

relação entre

b e T, segue que existe um valor para o parâmetro de impacto, b

min

, abaixo do

qual a transferência de energia será apenas

max

ao invés do valor dado pela equação

(2.14), ou seja,

min

= . ( 2.18 )

Para um elétron livre, não existe limite inferior para

T, mas para elétrons ligados em um

átomo, é necessário considerar seu movimento. Se o “tempo de colisão”, )/(

vb=

(isto

assume, razoavelmente, que a maior parte da interação entre a partícula e o elétron

acontece a uma distância b ao longo do eixo x), for comparável ao período orbital (

/1 )

do elétron, a partícula produzirá uma deformação na órbita deste elétron, mas não irá

transferir energia para ele. Assim, se for definido )/(

max

vvb

, nenhuma energia será

transferida para qualquer valor de

b maior do que b

max

. A relação entre b e T, equação

(2.14) indica que a menor transferência de energia possível é dada por:

242

min

== . ( 2.19 )

A seção de choque diferencial por elétron,

d , é definida como sendo a área em que a

partícula perderá uma quantidade de energia entre

T e dTT

, que é igual à área de um

anel de raio

b e largura db; então

dbbd 2

, ( 2.20 )

que em termos da energia transferida (

T) fica:

= . ( 2.21 )

Como as partículas incidentes são distribuídas com igual probabilidade em toda região em

torno do elétron, a probabilidade para uma colisão com transferência de energia entre

T e

dTT + será proporcional a

d e portanto irá decrescer muito rapidamente com o aumento

da energia transferida e irá também diminuir com a velocidade

v da partícula incidente.

Agora é, então, possível determinar a perda de energia média

∆

para uma partícula que

percorre uma distância

x∆ em um material com uma densidade de N átomos por unidade

de volume e com

Z elétrons por átomo (assumindo que

∆

). O número de colisões

dn em que ocorre uma perda de energia entre T e dTT

é igual ao número total de

elétrons por unidade de área,

, multiplicado pela seção de choque diferencial

d ;

ddn = . A perda de energia dentro deste intervalo de energia, devido a todas as

colisões, será então

dTdnT = . A perda de energia total

∆

ao longo de x∆ é obtida

integrando-se sobre todas perdas de energia diferenciais possíveis.

∫∫

==∆

max

min

dTdnTE

συ

( 2.22 )

onde

d depende da energia [equação. (2.21)]. Substituindo

d em (2.22):

⎟

⎠

⎞

⎜

⎝

⎛

==∆

∫

min

max

min

Tmv

υπυπ

. ( 2.23 )

Substituindo

xNZ ∆=

e as equações (2.17) e (2.19) e dividindo por x ∆ :

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

∆

242

NZez

. ( 2.24 )

Este resultado é válido para

zv pequeno. É possível também determinar uma expressão

para

xE ∆∆ / utilizando-se a distribuição de elétrons no material: A partícula incidente

perderá uma mesma quantidade de energia

T para todos elétrons que estiverem a uma

mesma distância

b da sua trajetória. O número dn de elétrons recebendo energia entre T e

dTT + inclui todos elétrons localizados dentro de um cilindro de diâmetro b, espessura db

e comprimento

∆

: xdbbNZdn

∆

= 2

. A quantidade total de energia absorvida por

todos estes elétrons é

(

)

(

)

xbdbmvNZezdnT ∆= 4

242

[equação. (2.14)]. Assim, a

energia total perdida para todos elétrons no material é obtida simplesmente integrando esta

expressão em

b (somando as contribuições de cada cilindro):

⎟

⎠

⎞

⎜

⎝

⎛

∆

min

max

NZez

. ( 2.25 )

O fator multiplicado pelo logaritmo é o mesmo na teoria quântica, mas os valores

escolhidos para

min

b e

max

b são diferentes. Os valores clássicos foram obtidos

anteriormente e a razão entre eles dá

)()(

minmax

zemvbb = , que substituído na equação

(2.25) leva ao mesmo resultado de antes, [equação. (2.24)].

Para partículas incidentes com velocidades mais altas, é necessário considerar a natureza

ondulatória das partículas, uma vez que o parâmetro de impacto não pode ser menor do que

o tamanho do pacote de onda que representa a partícula. Considerando a colisão no sistema

de coordenadas do centro de massa, a partícula pesada,

M, estará em repouso, e o elétron

estará movendo-se ao longo do eixo

x com velocidade v na direção de M. Após a colisão, o

elétron terá componentes do momento em ambas direções,

x e y. O maior valor possível

para o momento na direção

y será igual a mv (correspondente a uma deflexão de 90º). O

princípio da incerteza determina que

h~ ≥

∆

py . O que significa que a definição do

parâmetro de impacto dentro de um intervalo

∆

produzirá uma incerteza

∆

momento

p , onde

()

∆

≥∆ h~ . Como o momento na direção y não pode ser maior

que

mv, a incerteza máxima não pode exceder mv; conseqüentemente, o menor valor de

y∆ será aproximadamente

(

)

mvy h~

min

∆ , que no sistema de coordenadas do laboratório

fica

(

)

mvy 2~

min

h∆ . Isto, por sua vez, significa que o parâmetro de impacto não pode

ser menor do que

min

y∆ ; assim, o valor

min

= ( 2.26 )

é obtido no tratamento quântico em contraste com o valor clássico

(

)

min

mvzeb = . Para

substituir na equação (2.25) deve-se escolher o maior destes valores; assim, é conveniente

definir a razão,

(

)

vze h2

, entre os valores clássico e quântico para

min

b . Desta

forma, o resultado clássico para

∆

será válido para 1>>

, ao passo que para 1

deverá ser utilizado o resultado da mecânica quântica:

⎟

⎠

⎞

⎜

⎝

⎛

∆

4 mv

NZez

. ( 2.27 )

h é a energia de excitação do elétron, portanto é necessário o uso de um valor médio,

, abrangendo todos valores possíveis de

No tratamento quântico, Bethe mostrou que a seção de choque é a mesma da equação

(2.21) exceto pelo fato de

T ter que ser substituído por mQ 2/

g= , onde 'p-pg = é a

variação do momento da partícula incidente. Para um elétron livre,

g é igual ao q definido

na equação (2.11) já para um elétron ligado, uma maior parcela de

g pode ser transferida ao

átomo. Além disto, a equação (2.21) deve ser multiplicada pelo “fator de forma”

)(g

(Fano, 1963 citado por ATTIX, 1968), que representa a probabilidade de excitação de um

átomo para o

i-ésimo estado (incluindo estados contínuos) quando este absorve um

momento

g. Conseqüentemente,

min

Q deverá ser determinado como função da energia de

excitação

U do i-ésimo nível, ficando

min

2mvUQ

= ;

max

Q permanecerá o mesmo:

max

2mvQ = . Na integral da equação (2.23), )/ln(

minmax

TT é equivalente a

minmax

)/2()/ln(

UmvQQ = . Após a soma sobre todos os níveis U

, a energia de excitação

média

I como definida na equação (2.29) é obtida, e, incluindo termos relativísticos, o

resultado quântico completo fica:

⎭

⎬

⎫

⎩

⎨

⎧

−−

⎥

⎦

⎤

⎢

⎣

⎡

−

∑

NZez

)1(

, ( 2.28 )

com )/(

cv=

. Os fatores de forma realizam as correções relativas às camadas eletrônicas

e, portanto, não é necessária a inclusão destas correções na teoria.

A energia de excitação média do átomo

I é definida por

∑

UfIZ lnln , ( 2.29 )

onde

é a força do oscilador (oscillator strength) para a transição com energia de

excitação

. O termo

()

∑

representa as correções para camadas, que se fazem

necessárias devido ao fato de os elétrons contribuírem menos para o poder de frenamento

quando

v é comparável às suas velocidades orbitais. Para velocidades baixas, C

é bastante

grande, mas assume valores absolutamente pequenos para altas velocidades, sendo

inversamente proporcional a

em uma aproximação de primeira ordem.

Elétrons

Enquanto passa através de um material, elétrons experimentam interações similares às

interações de partículas pesadas, porém, duas diferenças básicas se manifestam: (a) nas

colisões com elétrons atômicos, podem ocorrer grandes perdas de energia; e (b) elétrons

com energias de poucas centenas de keV’s apresentam efeitos relativísticos.

As expressões para partículas pesadas são válidas para elétrons de baixa energia com a

simples modificação do parâmetro de impacto mínimo da equação (2.26). Como o poder

de frenamento não-relativístico é pouco utilizado quando se trata de elétrons, será

apresentada apenas a expressão relativística. Devido a indistinguibilidade dos elétrons,

assume-se que após uma colisão com outro elétron, o elétron primário é aquele com maior

energia. A perda de energia média para elétrons (e pósitrons) por colisão é então

()

/gMeV.cm ,

1535.0

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

−∆+

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

−

δτ

βρ

mcI

col

( 2.30 )

onde, para elétrons

()

2 1=∆ e

(

)

[

]

(

)

[

]

() ( ) ()

[]

()

1 ln12 2 1

ln1

∆−++∆

∆−+∆∆−+−−=

−

ττ

τττβ

( 2.31 )

e para pósitrons

∆

()

(

)

()

(

)

(

)

(

)

()

()()() () ()

()

⎥

⎦

⎤

∆+∆−∆++

−

⎢

⎣

⎡

∆−∆++

∆

−∆+−∆=

434

4 1 3 4 131

3 1 31

τττ

ττ

( 2.32 )

Aqui,

(

)

mcT=

é a correção relativa à densidade do meio, e eV 006,511

=mc ;

∆

é a energia máxima entregue a raios δ dividido por

mc . Escrever estas equações contendo

∆ explicitamente é muito útil quando se deseja uma descrição em função do LET

(ATTIX, 1968).

2.2.2 Alcance

Com as constantes colisões e eventual emissão de radiação de bremsstrahlung, as

partículas carregadas penetram num meio material até que sua energia cinética entre em

equilíbrio térmico com as partículas do meio, estabelecendo um “alcance” no meio

absorvedor, após um percurso em linha reta ou em zig-zag. As partículas pesadas, como

alfa e fragmentos de fissão, têm uma trajetória praticamente retilínea dentro do material, ao

contrário dos elétrons que têm uma trajetória tortuosa. Para cada tipo de partícula pode-se

definir um alcance, utilizando variações da definição, devido às dificuldades experimentais

em sua determinação (TAUHATA, 2002).

Uma partícula carregada percorre um caminho aproximadamente bem definido no material

antes de dissipar toda sua energia. A distância ao longo deste caminho, ou seja, a soma de

todos os deslocamentos desde o ponto inicial até o final, é chamada de

comprimento da

trajetória

da partícula. Devido à natureza estatística do processo de perda de energia,

partículas com uma mesma energia inicial

E terão trajetórias de comprimentos

ligeiramente diferentes (fenômeno denominado dispersão do alcance ou “range

straggling”). A média destes vários comprimentos de trajetória é definida como

alcance

. Os comprimentos das trajetórias das partículas que sofreram interações catastróficas,

como ocorre nas reações nucleares ou na emissão bremsstrahlung, não são considerados

nesta média.

Muitas vezes, não é possível observar a trajetória das partículas, mas apenas a projeção

desta em uma dada direção. Isto ocorre, por exemplo, quando um feixe de partículas incide

perpendicularmente em placas de diferentes espessuras de um dado material e observa-se o

número de partículas que o atravessam. Se )(

tN é o número de partículas que atravessam

uma placa de largura

t, o alcance projetado é definido por

∫

∞

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

)(

tdN

( 2.33 )

onde

N é o número de partículas incidentes na placa excluindo aquelas que sofrem

eventos catastróficos. O alcance projetado é então a distância média que partículas de uma

mesma energia inicial penetram em um material na direção de incidência. Outras

grandezas úteis são o

alcance projetado médio (

médio

t ) – valor de t em que 21)(

tN , ou

seja, é o ponto em que metade do número original de partículas é absorvido/transmitido, e

alcance projetado mais provável ( 't ) – valor de t para o qual dttdN )( é máximo.

A FIGURA 2.4 mostra uma curva típica para o alcance de partículas

em um gás. Esta

curva pode ser obtida usando-se uma separação fixa entre a fonte de partículas

e o

detector e mudando a pressão do gás entre eles. Nela está plotado o número

N de partículas

detectadas em um certo intervalo de tempo versus a espessura,

t, atravessada no gás, que

certamente é proporcional à pressão deste. A região de dispersão ocorre entre

aproximadamente

t95.0 e t05.1 , sendo que muito poucas partículas são perdidas entre

0 e t95.0 . Estas perdas são causadas principalmente por reações nucleares e

espalhamento múltiplo. Para partículas altamente energéticas, as reações nucleares

eliminam uma grande fração das partículas incidentes. Na FIGURA 2.4estão representados

o alcance projetado e o alcance projetado médio.

Na FIGURA 2.5 está representada uma curva de alcance para elétrons obtida com o

mesmo arranjo experimental. Neste caso, os freqüentes e grandes espalhamentos levam a

uma distribuição com profundidades de penetração muito diferentes. Por isso, o conceito

de alcance projetado não é muito útil para elétrons porque eles se espalham tanto que as

curvas )(

tN não possuem uma região de decrescimento rápido como nas curvas para

Partículas transmitidas

Espessura do absorvedor

médio

FIGURA 2.4: Alcance projetado para partículas α monoenergéticas. A assimetria

após o ‘joelho’ da curva é exagerada para mostrar mais claramente a diferença entre o

alcance projetado

t , e o alcance projetado médio,

médio

t (ATTIX, 1968).

Partículas transmitidas

FIGURA 2.5: Absorção de elétrons monoenergéticos. R

é definido como o

alcance extrapolado (ATTIX, 1968).

partículas carregadas pesadas. Sendo mais difícil a definição do alcance projetado, define-

se o

alcance extrapolado,

R , como indicado na FIGURA 2.5.

Nem o ponto final da trajetória da partícula, nem a espessura de material que esta é capaz

de penetrar são grandezas bem definidas. Por exemplo, uma partícula α, após termalizar, se

tornará, presumidamente, um átomo de hélio neutro que certamente irá se difundir através

do material indefinidamente. Portanto, é necessário definir um limite inferior para a

energia da partícula no qual considera-se que esta alcançou o ponto final da trajetória. Este

pode ser, por exemplo, a energia na qual a partícula não é mais capaz de produzir efeitos

em tecidos biológicos, ou ainda o limiar de energia em que o detector é capaz de identificar

a presença da partícula. Outro fator que prejudica a definição do alcance projetado é o fato

de que, uma vez que as partículas emergem do material absorvedor em várias direções,

nem todas elas são detectadas. Esta distribuição da direção em que a partícula deixa o

material resulta dos pequenos espalhamentos Colombianos. Estes espalhamentos são

também a origem das trajetórias tortuosas dentro do material, que por sua vez são

responsáveis pela diferença entre o alcance

R e o alcance projetado t . Este efeito é

conhecido como espalhamento múltiplo. A diferença

tR− geralmente é pequena para

partículas pesadas. Estes conceitos de comprimento da trajetória, alcance e alcance

projetado deixam de ser úteis em altas energias onde um número considerável de partículas

começa a sofrer interações catastróficas ao invés de sofrer perdas de energia

aproximadamente contínuas, através de excitação e ionização. Em energias muito altas é

mais útil tratar o movimento de prótons através da matéria da forma como é feita para

fótons e nêutrons, isto é, utilizando coeficientes de atenuação. Da mesma forma, elétrons

de altas energias geram cascatas de fótons e elétrons que dificultam bastante a definição do

conceito de alcance.

Uma grandeza teórica muito útil é o alcance CSDA (“continuous-slowing-down-

approximation”) que assume que a partícula perde energia continuamente ao longo de sua

trajetória e é definido por

)(

ERdE

∫

( 2.34 )

para todas as partículas carregadas. Como a expressão de Bethe para o poder de

frenamento não é válida para pequenas energias, a integração para

R começa em uma

certa energia

E e um alcance CSDA residual, experimentalmente estimado, é adicionado

à integral. O

alcance CSDA não é igual ao alcance

porque a atenuação das partículas é

discretizada em pequenas (mas não infinitesimais) perdas de energia. A diferença entre

eles é de 0.2% ou menos para prótons; para elétrons espera-se que esta diferença seja

maior.

alcance máximo,

máx

R , é definido como a espessura de material em que todo o feixe de

elétrons, inicialmente monoenergético e paralelo, é completamente absorvido

(comprimento da trajetória máximo). Este alcance é difícil de ser determinado para

elétrons com

≤

E keV, enquanto que para elétrons com 40 keV ≤ E ≤ 160 keV, ele pode

ser reproduzido mais precisamente do que os alcances definidos por outros métodos

(ATTIX, 1968).

3 Técnicas de Monte Carlo para o transporte de elétrons e fótons

O transporte de radiação na matéria tem sido objeto de intensos trabalhos desde o começo

do século XX. Hoje, sabe-se que fótons, elétrons e pósitrons, ao penetrarem a matéria,

sofrem múltiplas interações através das quais energia é transferida para os átomos e

moléculas do material e partículas secundárias são produzidas. Por meio de repetidas

interações com o meio, uma partícula origina uma cascata de partículas que é usualmente

chamada de chuveiro. Após cada interação, a energia da partícula é reduzida e novas

partículas são geradas de forma que a evolução de um chuveiro representa uma degradação

efetiva da energia.

Na técnica de Monte Carlo para a simulação do transporte de radiação, a história

(percurso) de uma partícula é uma seqüência aleatória de deslocamentos que terminam

com uma interação onde a partícula muda sua direção de movimento, perde energia e,

ocasionalmente, produz partículas secundárias. Para simular estas histórias é necessário um

modelo de interação, isto é, um conjunto de seções de choque diferenciais (SCD) dos

mecanismos de interação relevantes. As SCD’s determinam as funções distribuição de

probabilidade (FDP) das variáveis aleatórias que caracterizam a trajetória: 1) livre caminho

entre dois eventos de interação sucessivos, 2) tipo de interação, 3) perda de energia e

deflexão angular em um dado evento (e o estado inicial das partículas secundárias

produzidas, se for o caso). Uma vez conhecidas as FDP’s, as histórias aleatórias podem ser

geradas utilizando-se métodos de amostragem apropriados [incluindo a geração de

números (pseudo-)aleatórios]. Se o número de histórias é grande o suficiente, pode-se obter

informações quantitativas do processo de transporte simplesmente fazendo uma média de

valores obtidos nas histórias simuladas (SALVAT, 2001). Uma vez que a história de uma

partícula pode ser acompanhada individualmente, esta técnica pode ser usada para obter

informações a respeito de flutuações estatísticas de tipos de eventos específicos, além

disto, pode-se também obter respostas impossíveis de se obter por meio de experimentos,

como por exemplo: “Qual fração destes elétrons foi gerada no colimador e no filtro?” ou

“Qual a fração dos fótons que sofreu espalhamento Compton?” (ROGERS, 1990).

Inicialmente, o estudo de problemas de transporte de radiação era feito com base na

equação de transporte de Boltzmann. Porém, este procedimento apresentava consideráveis

dificuldades quando aplicado a geometrias limitadas – os métodos numéricos baseados na

equação de transporte tinham certo sucesso somente para geometrias simples,

essencialmente para meio infinito e semi-infinito (ZHENG-MING E BRAHME, 1993

citado por ROGERS, 1990). No final da década de 50, com a disponibilidade de

computadores, os métodos de simulação de Monte Carlo foram desenvolvidos como uma

poderosa alternativa para lidar com problemas de transporte (ROGERS, 1990). O método

Monte Carlo fornece a mesma informação que a solução da equação de transporte de

Boltzmann com o mesmo modelo de interação (seções de choque), mas é mais fácil de

implementar (BERGER, 1963 citado por SALVAT, 2001) e é capaz de responder a

questões relacionadas à história de uma partícula específica de maneira mais direta

(ROGERS, 1990). A principal desvantagem do método Monte Carlo está relacionada à sua

natureza aleatória. Todos os resultados são afetados por incertezas estatísticas, que podem

ser reduzidas aumentando-se o número de amostras e, conseqüentemente, aumentando o

tempo computacional. Em certas circunstâncias, as incertezas estatísticas podem ser

reduzidas utilizando-se técnicas de redução de variância (RUBISTEIN, 1981; BIELAJEW,

ROGERS, 1988 citado por SALVAT, 2001). Muitos problemas em dosimetria da

radiação, física de radioterapia e proteção radiológica vêm sendo tratados através de

técnicas de Monte Carlo porque a complexidade do transporte de elétrons e fótons no meio

material torna as soluções analíticas intratáveis (ROGERS, 1990).

A simulação do transporte de elétrons é muito mais difícil do que a de fótons. A principal

razão disto é que a energia média que um elétron perde em cada interação é muito pequena

(da ordem de poucas dezenas de eV). Conseqüentemente, elétrons de alta energia sofrem

um grande número de interações até serem efetivamente absorvidos pelo meio. Na prática,

simulações detalhadas são factíveis apenas quando o número médio de colisões não é

muito grande (até poucas centenas). Situações experimentais tratáveis através de simulação

detalhada são aquelas que envolvem ou fontes de elétrons com energia cinética inicial

baixa (até aproximadamente 100 keV) ou geometrias especiais como um feixe de elétrons

incidindo em uma lâmina fina. Para energias iniciais grandes e geometrias maiores e/ou

mais complexas, o número médio de colisões que um elétron sofre até ser absorvido torna-

se muito grande e a simulação detalhada ineficiente (SALVAT, 2001). Em função disto,

faz-se uso da técnica de história condensada (BERGER, 1963 citado por ROGERS, 1990),

em que o trajeto do elétron é dividido em uma série de passos no qual os efeitos de um

grande número de interações individuais são agrupados. Um grupo é composto pelas várias

deflexões decorrentes do espalhamento elástico; este agrupamento é feito utilizando-se

uma teoria de espalhamento múltiplo como a de Molière (1948) ou de Goudsmit e

Saunderson (1940) (citado por ROGERS, 1990). O outro usa um modelo de ‘frenamento

contínuo’ para agrupar o grande número de pequenas perdas de energia.

Os códigos de Monte Carlo podem ser divididos em duas categorias, chamadas de classe I

e classe II por Berger (1963) (citado por ROGERS, 1990). Eles se distinguem pela maneira

como tratam os eventos que levam à produção de fótons de bremsstrahlung e/ou a

produção de elétrons secundários. Nos modelos de classe I, as perdas de energia e

deflexões angulares associados a todos eventos são agrupados e a energia e direção do

elétron primário não são afetados pela criação de partículas secundárias. Nos modelos de

classe II, as interações afetam a energia e a direção de movimento do elétron primário, mas

apenas quando os elétrons secundários ou fótons criados têm energia acima de um certo

valor limite; sendo ainda agrupados os efeitos da produção de partículas secundárias com

energia abaixo deste limite.

Um código de simulação Monte Carlo possui quatro componentes principais: (1) dados de

seções de choque para todos os processos considerados na simulação, (2) algoritmos

usados para o transporte das partículas, (3) métodos usados para especificar a geometria do

problema e para calcular as grandezas físicas de interesse e (4) análise das informações

obtida durante a simulação. Os dois primeiros podem afetar consideravelmente o tempo

computacional, porém eles não afetam o aspecto físico da simulação (ROGERS, 1990).

3.1 Modelos classe I e classe II

Berger (1963) (citado por ROGERS, 1990) dividiu os algoritmos de transporte de elétrons

em duas principais classes diferenciadas pela forma com que a energia do elétron primário

se relaciona com a energia perdida nas interações com o meio. Nos modelos de classe I, em

cada passo da história condensada os efeitos de todas interações de um certo tipo são

agrupados. Modelos de classe II agrupam apenas os efeitos de um subconjunto de

interações de um dado tipo e trata os efeitos das interações restantes individualmente,

motivo pelo qual é chamado também de simulação mista. Por exemplo, para perdas de

energia por colisão, é usado um modelo de perda de energia contínua que agrupa os efeitos

de todas interações que produzem elétrons com energia abaixo de algum limite arbitrário.

Trata-se individualmente apenas as relativamente raras interações catastróficas que criam

partículas secundárias acima deste mesmo limite arbitrário de energia. Essas interações

fazem com que o elétron primário perca energia e seja defletido. A escolha dos limites de

energia para se considerar a criação de elétrons secundários ou fótons como eventos

discretos é arbitrária e é um componente do algoritmo, e não do processo físico envolvido.

Uma forma simples de cálculo de Monte Carlo para elétrons é o modelo de ‘aproximação

do frenamento contínuo’ (continuous-slowing-down-approximation – CSDA) no qual não

se considera a produção de partículas secundárias e utiliza-se o poder de frenamento total

irrestrito (apêndice B) para o cálculo da perda de energia em cada passo. Como o próprio

nome diz, considera-se que o elétron perde energia continuamente ao longo de sua

trajetória, apesar do fato desta trajetória ser dividida em um número finito de passos

fazendo parecer que a energia é degradada de maneira discretizada. O modelo CSDA é

claramente um algoritmo de classe I. Contudo, algoritmos de classe I podem ser

sofisticados o suficiente para incluir a geração de partículas secundárias e levar em conta a

dispersão da perda de energia (energy-loss straggling). A parte da perda de energia

PASSO SEGUIDO DE ESPALHAMENTO M

LTIPLO

EVENTO

DISCRETO. O

ELÉTRON CRIADO

É SEGUIDO

SEPARADAMENTE

A ENERGIA

PERDIDA CONTINUAMENTE

PARA ELÉTRONS E RAIOS-X DE BAIXA

ENERGIA

FIGURA3.1: Representação das perdas de energia e deflexões de um elétron em um algoritmo de classe II.

Considera-se que a energia é depositada continuamente ao longo dos passos, apesar do fato de as partículas

secundárias produzidas com energias menores do que o valor de corte difundirem esta energia em uma região

em torno da trajetória (área sombreada). No final de cada passo o elétron é defletido de um pequeno ângulo e

função do espalhamento múltiplo. Nos eventos de interação discretos são criados elétrons ou fótons de

bremsstrahlung com energias acima do valor de corte.

contínua é modelada usando o chamado poder de frenamento restrito (apêndice B). Esta é a

parte do poder de frenamento restrita à criação de partículas secundárias com energias

menores do que o limite considerado para produção de partículas secundárias.

A FIGURA 3.1representa o transporte de elétrons em um algoritmo classe II. Neste

modelo, o elétron se move em pequenos passos retilíneos. Em cada passo, é usada uma

teoria de espalhamento múltiplo para determinar o ângulo em que o elétron será defletido.

A teoria de espalhamento múltiplo deve também ser utilizada na determinação do

comprimento real da trajetória do elétron durante o passo. Deste comprimento pode-se

deduzir a quantidade de energia perdida continuamente, incluindo os elétrons e fótons

criados com energia abaixo do limite pré-determinado. Embora a energia destas partículas

de baixa energia seja depositada por toda área sombreada em torno da trajetória, no modelo

considera-se que esta energia é depositada no próprio trajeto. Depois de percorridas as

distâncias determinadas pelas seções de choque apropriadas, o elétron interage e produz

uma partícula secundária com energia acima do limite de produção. A partícula secundária

é “armazenada” e a energia do elétron primário é reduzida de uma quantidade

correspondente à energia desta partícula secundária. Esta aproximação, em geral muito

boa, ignora a energia de ligação dos elétrons alvo e qualquer energia absorvida pelo núcleo

ou elétron quando um fóton de bremsstrahlung é criado (ROGERS, 1990).

3.2 Simulação do transporte de radiação

3.2.1 Modelo de espalhamento e funções distribuição de probabilidade

Considere uma partícula com energia

(energia cinética, no caso de elétrons e pósitrons)

movendo em um certo meio (gasoso, líquido ou sólido amorfo, onde as moléculas são

distribuídas aleatoriamente com densidade uniforme). A composição deste meio é

especificada pela sua fórmula estequiométrica, isto é, número atômico

Z e número de

átomos por molécula

n de todos elementos presentes. Os índices estequiométricos

n não

precisam necessariamente ter valores inteiros, e no caso de ligas, por exemplo, eles devem

ser iguais à porcentagem de cada elemento. O peso molecular é

∑

iiM

AnA , onde

A é o

peso atômico do

i-ésimo elemento. O número de moléculas por unidade de volume é dado

por

=Ν , ( 3.1 )

onde

N é o número de Avogadro e

a densidade do material.

Em cada interação, a partícula pode perder uma energia

W e/ou mudar sua direção de

movimento. A deflexão angular é determinada pelo ângulo polar de espalhamento

–

ângulo entre a direção da partícula antes e depois da interação – e pelo ângulo azimutal

Assumindo que a partícula pode interagir com o meio através de dois mecanismos

independentes “A” e “B” (por exemplo, espalhamento elástico e inelástico, no caso de

elétrons de baixa energia). O modelo de espalhamento é completamente especificado pelas

seções de choque diferenciais moleculares (SCD)

dWd

),;(e),;(

ΩΩ

, ( 3.2 )

onde

Ωd é um elemento de ângulo sólido na direção ),(

, e a dependência paramétrica

da SCD com a energia da partícula

foi feita explicitamente. Considerando que as

moléculas do meio são orientadas aleatoriamente, a SCD é independente do ângulo de

espalhamento azimutal, ou seja, a distribuição angular das partículas espalhadas é simétrica

em torno da direção de incidência. As seções de choque totais (por molécula) são

∫∫

Ω

θθπσ

) ,;( ) ( 2 )( WE

dWd

dsendWE

. ( 3.3 )

As FDP’s para perda de energia e para o ângulo polar de espalhamento em eventos de

espalhamento individuais são

) ,;(

)(

) ,;(

θπ

dWd

sen

WEp

Ω

. ( 3.4 )

Note que

ddWWEp

) ,;( dá a probabilidade (normalizada) de, em um evento de

espalhamento do tipo A, a partícula perder energia no intervalo ),(

dWWW + e ser

defletida na direção com ângulo polar no intervalo ),(

. O ângulo de espalhamento

azimutal em cada colisão é uniformemente distribuído no intervalo )2,0(

, isto é

)( =

p . ( 3.5 )

A seção de choque total para interação é

)()()( EEE

BAT

. ( 3.6 )

Quando a partícula interage com o meio, o tipo de interação que ocorre é uma variável

aleatória, que assume os valores “A” e “B” com as probabilidades

p p

== e

. ( 3.7 )

Vale lembrar que este tipo de modelo de espalhamento individual somente é válido quando

os efeitos da difração devido ao espalhamento coerente por vários centros são desprezíveis.

O que significa que a simulação é aplicável apenas a meios amorfos e, com algum cuidado,

a sólidos policristalinos.

Para um entendimento intuitivo do processo de espalhamento, pode-se imaginar cada

molécula como sendo uma esfera de raio

r de forma que a área de seção de choque

seja igual à seção de choque total

. Assumindo que uma partícula incide

normalmente em uma lâmina muito fina de largura

ds , o que a partícula incidente vê é

uma distribuição uniforme de

esferas por unidade de área; quando esta partícula

atinge uma destas esferas ocorre uma interação. Assim, a probabilidade de ocorrer uma

interação nesta lâmina de material é igual à fração da área ocupada pelas esferas,

Em outras palavras,

é a probabilidade de interação por unidade de distância

percorrida no meio. Sua inversa,

)(

−

Ν≡

σλ

, ( 3.8 )

é o livre caminho médio total entre interações.

Considerando agora uma partícula que move em um meio infinito. A FDP )(

sp para

comprimento da trajetória da partícula entre dois eventos de interação pode ser obtida

como a seguir: A probabilidade da partícula percorrer uma distância

s sem interagir é

∫

∞

=ℑ

' )'()(

dssps . ( 3.9 )

A probabilidade

dssp )( da próxima interação ocorrer quando a distância percorrida estiver

no intervalo ),(

dsss + é igual ao produto de )(s

ℑ

(probabilidade de chegar em s sem

interagir) e

1−

(probabilidade de interagir em ds ), isto é

∫

∞

−

' )'()(

dsspsp

. ( 3.10 )

A solução desta equação integral, com a condição de contorno 0)( =∞

p , é a familiar

distribuição exponencial

)exp()(

ssp

λλ

−=

−

. ( 3.11 )

Note que o livre caminho médio

coincide com o comprimento médio da trajetória entre

colisões sucessivas:

∫

∞

)(

dsspss

. ( 3.12 )

A diferencial do inverso do livre caminho médio para o processo de interação A é definido

como

),;(),;(

212

dWd

Ω

Ν=

Ω

−

. ( 3.13 )

Evidentemente, a integral desta diferencial dá o inverso do livre caminho médio para o

processo,

dWd

dsendW

σθ

θθπλ

Ν=

Ω

∫∫

−

) ,;( ) ( 2

. ( 3.14 )

O produto

é chamado de seção de choque macroscópica, embora este não seja um

nome apropriado para uma grandeza que tem dimensão de inverso de comprimento. Note

que o inverso do livre caminho médio total é igual à soma dos inversos dos livres caminhos

médios dos diferentes mecanismos de interação considerados,

111 −−−

BAT

λλλ

. ( 3.15 )

(SALVAT, 2001).

3.2.2 Geração das trajetórias aleatórias

A trajetória de cada partícula começa em uma dada posição, com direção inicial e energia

determinadas de acordo com as características da fonte. O estado de uma partícula

imediatamente depois de uma interação (ou após penetrar na amostra ou no começo da

trajetória) é definido pelas suas coordenadas ),,(

zyx

r , energia

e pelas componentes

do vetor unitário ),,(

wvu=d que representa os cossenos diretores da direção de

movimento. As coordenadas de posição e direção de movimento são definidas em relação a

um sistema de coordenadas retangular fixo, o referencial do laboratório que pode ser

definido arbitrariamente. Cada trajetória simulada é então caracterizada por uma série de

estados

r ,

E e

, onde

r é a posição em que ocorre o n-ésimo evento de

espalhamento e

E e

são a energia e os cossenos diretores da direção de movimento

exatamente após este evento. O comprimento

s do livre caminho até o próximo evento de

interação, o mecanismo de espalhamento envolvido, a mudança na direção e a perda de

energia nesta colisão são variáveis aleatórias que são “sorteadas” a partir das FDP’s

correspondentes usando os métodos descritos na seção 1.2. A geração das trajetórias

aleatórias é feita assumindo-se que trajetória já foi simulada até o estado

r ,

E ,

O comprimento do trajeto entre duas interações sucessivas é distribuído de acordo com a

FDP dada pela equação (3.11) Os valores aleatórios de

s são gerados a partir da fórmula

de amostragem

−

, ( 3.16 )

onde

representa um número aleatório uniformemente distribuído no intervalo (0,1). A

interação seguinte ocorre na posição

nnn

sdrr

. ( 3.17 )

O tipo desta interação (“A” ou “B”) é selecionado a partir das probabilidades dadas pela

equação (3.7) usando-se o método da transformada inversa (seção1.2.2). A perda de

energia

W e o ângulo de espalhamento polar

são obtidos através da distribuição

),;(

WEp

, equação (3.4) usando-se a técnica de amostragem mais conveniente. O

ângulo de espalhamento azimutal é gerado, de acordo com a distribuição uniforme em

(0,2π) [equação (3.5)], através da fórmula de amostragem

πξ

Após a determinação dos valores de

W ,

, a energia da partícula é reduzida,

WEE

−=

, e a direção de movimento após a interação )',','(

wvu

d é obtida

realizando a rotação de ),,(

wvu

=d (FIGURA 3.2). A matriz rotação ),(

R é

determinada pelos ângulos de espalhamento polar e azimutal. Para se obter o vetor

R dd

),(

φθ

após a interação deve-se observar que, se a direção inicial está ao longo

do eixo

z, a direção após a colisão é

FIGURA 3.2: Deflexões angulares em eventos de espalhamento.

)()(

cos

θφ

φθ

RRsensen

sen

⎟

⎠

⎞

⎜

⎝

⎛

, ( 3.18 )

onde )1,0,0(

z e

⎟

⎠

⎞

⎜

⎝

⎛

−

θθ

cos0sen

010

sen0cos

)(

⎟

⎠

⎞

⎜

⎝

⎛

−

100

0cossen

0sencos

)(

φφ

( 3.19 )

são as matrizes rotação que correspondem à realização das rotações em torno dos eixos

y e

z de ângulos iguais a

, respectivamente. Se, por outro lado,

forem,

respectivamente, os ângulos polar e azimutal da direção inicial

)cos,sensen,cos(sen

ϑϕϑϕϑ

d , ( 3.20 )

a rotação )()(

−−

RR transforma o vetor

em z

. Fica claro então que a direção final

d pode ser obtida realizando-se a seguinte seqüência de rotações no vetor inicial: 1)

)()(

−−

RR , que transforma

em z

; 2) )()(

RR , que gira o vetor z

de acordo

com os ângulos de espalhamento polar e azimutal ‘sorteados’; e 3) )()(

RR , que

inverte a rotação realizada no primeiro passo. Assim,

)( )( )( )( )( )(),(

−

zyyzyz

RRRRRRR . ( 3.21 )

A direção final é então

⎟

⎠

⎞

⎜

⎝

⎛

φθ

ϑϕφθ

cos

)()(

),(

sensen

sen

RRR

yznn

( 3.22 )

e seus cossenos diretores:

[]

. cos1cos'

, cos

cos'

, cos

cos'

φθθ

φφ

senwww

usenvw

sen

vsenuw

sen

−−=

−

( 3.23 )

Estas equações são indeterminadas quando

≈

w , ou seja, quando a direção inicial é

aproximadamente paralela ou anti-paralela ao eixo

z; nestes casos fica estabelecido que

cos,sensen,cossen

±= wv u . ( 3.24 )

Além disto, as equações (3.23) não são muito estáveis numericamente e a normalização de

d tende a ser diferente de 1 após repetido uso delas. Isto é remediado através da

renormalização periódica de

d .

A simulação da trajetória prossegue então repetindo-se estes passos e termina quando a

partícula deixa o material ou quando sua energia torna-se menor do que a energia de

absorção

abs

E , na qual assume-se que a partícula foi efetivamente freada e absorvida pelo

meio (SALVAT, 2001).

3.2.3 Seleção de um processo físico

Praticamente todos os processos envolvidos no transporte de uma partícula têm natureza

aleatória; na melhor das hipóteses, pode-se conhecer as distribuições de probabilidade que

governam um evento. Por exemplo, o conhecimento da seção de choque total de um fóton

em um material não implica em saber o quanto o fóton vai penetrar em um meio, mas

apenas o valor médio do comprimento de sua trajetória até a próxima interação. Da mesma

forma, seções de choque diferenciais fornecem a apenas a probabilidade de um processo

ocorrer como função de alguma variável representativa do estado final da partícula (por

exemplo, energia ou ângulo de espalhamento). Portanto, um elemento essencial em

qualquer simulação de Monte Carlo é a capacidade de realizar a amostragem a partir de

várias distribuições de probabilidade que descrevem a física dos processos envolvidos e a

simulação da natureza aleatória destes eventos individuais. Esta pode ser uma tarefa muito

difícil, mas, felizmente, já existem algoritmos precisos e eficientes para amostragem a

partir das distribuições mais utilizadas. Para os processos menos freqüentes, introduz-se

aproximações que tornam as rotinas de amostragem ineficientes, mas como elas são usadas

com pouca freqüência, seu efeito no tempo e na exatidão da simulação é pequeno.

Veja um exemplo. Em um modelo de transporte de fótons que considera apenas

espalhamento Compton e eventos de produção de pares, a seção de choque total em uma

dada energia é dada por

ppComptontot

( 3.25 )

As perguntas em questão são: qual distância um fóton caminha antes de interagir e qual

interação ocorre? O procedimento começa com a escolha de dois números aleatórios

uniformemente distribuídos entre 0 e 1. Sabendo que os comprimentos das trajetórias

dos fótons são exponencialmente distribuídos, faz-se uso de uma distribuição exponencial

com um caminho médio dado por

1−

tot

para amostragem. De acordo com as informações

acima, a variável

(cm) x

tot

−= ( 3.26 )

é exponencialmente distribuída entre zero e infinito com valor médio

tot

Σ1 . Desta forma,

determina-se que um fóton percorrerá uma distância igual a

cm para então interagir.

Pode-se então determinar que tipo de interação ocorre escolhendo uma interação Compton

totCompton

Σ≤

e uma interação de produção de pares caso contrário. Estas são as

duas rotinas mais simples e mais freqüentemente utilizadas em simulações de Monte Carlo

(ROGERS, 1990).

3.2.4 Efeitos do tamanho do passo da partícula

A escolha do tamanho do passo em simulações do transporte de elétrons/pósitrons pode

afetar dramaticamente a precisão dos resultados e o tempo computacional.

No princípio do passo de um elétron, sabe-se a posição,

r , e a direção,

, iniciais, sendo

necessária, então, a escolha do comprimento total do caminho

s até o final do passo. Este

comprimento total do caminho é usado por várias teorias físicas para determinar

(possivelmente estocasticamente) a energia do elétron no final do passo, a posição

r e a

nova direção

no ponto final, como caracterizados pelo ângulo de espalhamento

múltiplo

Para calcular a posição do ponto final

r do passo, é necessário calcular s, a componente

do deslocamento ao longo da direção inicial

. As grandezas t e s estão representadas na

FIGURA 3.3 para um passo típico de um elétron. A grandeza

s se relaciona com

r através

da expressão

FIGURA 3.3: Componentes geométricas do passo de um elétron.

iif

s drr

)( −= . ( 3.27 )

Berger propôs a seguinte relação entre

t e s :

)]( cos1[

tts

+= . ( 3.28 )

Note que

s, como calculado pela equação (3.30), está correlacionado com o ângulo de

espalhamento múltiplo e portanto exibe uma distribuição em torno do seu valor médio.

Lewis (1950) (citado por ROGERS, 1990) determinou a expressão exata para

s ,

∫

tdts

)'( cos'

. ( 3.29 )

A “correção para comprimento do caminho” (CCC), definido como

(

)

sst − , é uma

função dependente de

t e é uma medida da curvatura no passo do elétron. Para ilustrar isto,

a FIGURA 3.4 apresenta a CCC em água versus a energia cinética do elétron para vários

tamanhos de passo, caracterizados pela fração da energia cinética do elétron perdida em

todos processos de colisão. A partir desta figura pode-se concluir que a CCC é significante

exceto para elétrons de energia muito alta ou para passos muito pequenos. Assim, para

cálculos precisos deve-se incluir a correção para a curvatura do caminho, a menos que se

aceite pagar o preço de calcular histórias usando passos extremamente curtos.

A FIGURA 3.4 ilustra alguns fatos interessantes a respeito do transporte de elétrons. À

medida que a energia do elétron aumenta, os elétrons tendem a deslocar em linha reta, e a

FIGURA 3.4: Correção para o comprimento do passo versus a energia

cinética para vários tamanhos de passo caracterizados pela perda percentual de

energia devido a todos eventos de colisão.

curvatura do caminho decorrente do espalhamento múltiplo torna-se menos importante. À

medida que a energia do elétron decresce, nota-se um nivelamento das curvas com perda

de energia constante por passo. Isto reflete o fato de que, embora os passos do elétron

tenha mais curvatura em baixas energias, o poder de frenamento de colisão também

aumenta, reduzindo o tamanho relativo do passo para uma dada perda de energia.

É necessário considerar também o deslocamento lateral

do elétron durante o percurso

do passo como ilustrado na FIGURA 3.3. Berger (1963) propôs a expressão

θρ

sen

t= ( 3.30 )

que expressa a correlação básica de

com o ângulo de espalhamento múltiplo

. Assim

como no caso da CCC, é necessário incluir o deslocamento lateral (DL) para uma

descrição precisa das trajetórias no transporte de elétrons, exceto em altas energias ou para

passos muito pequenos.

Os resultados de cálculos que consideram a CCC e o DL independem do tamanho do passo

se o transporte do elétron estiver sendo feito com precisão. O efeito do deslocamento

lateral é menos significante e, para passos de comprimento pequeno, os resultados que não

consideram a CCC e o DL convergem para os resultados com CCC e DL, refletindo o fato

de que passos de menor tamanho determinam um transporte mais preciso. Ao invés de

considerar as correções da curvatura e deslocamento lateral do passo do elétron, pode-se

obter resultados precisos reduzindo o tamanho passo do elétron. Porém, deve-se tomar

cuidado para não reduzir demais o passo violando os limites da teoria de espalhamento

múltiplo; caso contrário, podem ocorrer artefatos computacionais. É mais desejável incluir

a CCC e o deslocamento lateral, através dos quais a dependência com o tamanho do passo

é bastante reduzida. Isto permite o uso de passos maiores e tempos de execução

potencialmente mais rápidos.

É necessário um maior cuidado ao simular o transporte de elétrons na vizinhança de

interfaces. A razão disto é que as teorias de espalhamento múltiplo usadas em simulações

de Monte Carlo de história condensada são válidas somente em meio infinito ou semi-

infinito. A introdução de interfaces viola limitações fundamentais destas teorias. Para

evitar os artefatos dos cálculos faz-se necessário, pelo menos na vizinhança de uma

interface, encurtar os passos do elétron de maneira que para a maioria deles o transporte

ocorra em um meio infinito. Qualquer que seja a estratégia empregada, a razão

fundamental é a mesma: a maioria dos trajetos dos elétrons deve acontecer em um meio

infinito. Isto é consumado diminuindo-se o tamanho dos passos, pelo menos na vizinhança

de interfaces. Quando um elétron atravessa uma interface ocorre alguma violação das

teorias envolvendo espalhamento múltiplo, mas como isto afeta apenas uma pequena

fração da trajetória total do elétron a qualidade dos resultados não é comprometida

desnecessariamente (ROGERS, 1990).

3.2.5 Transporte de partículas como um processo Markoviano

Num processo Markoviano os valores que uma variável aleatória assumirá ao longo da

simulação são determinados estatisticamente através de uma relação com eventos atuais e

dependem apenas do evento imediatamente anterior. Devido à característica Markoviana

do transporte de partículas, a história de uma partícula pode ser interrompida em um estado

arbitrário (qualquer pondo da trajetória) e ser retomada a partir deste estado sem qualquer

influência nos resultados.

Em simulações mistas (classe II) de transporte de elétrons/pósitrons, é necessário limitar o

tamanho

s de cada passo da partícula de forma que este não exceda um dado valor

max

s .

Para se fazer isto, obtém-se o livre caminho

s até a próxima interação a partir da mesma

FDP exponencial (3.11), mas quando

max

ss > , à partícula é permito deslocar apenas a

distância

max

s ao longo da direção de movimento. No final deste passo truncado nada é

feito (a partícula mantém sua energia e direção de movimento inalteradas); porém, por

conveniência, assume-se que a partícula sofreu uma interação delta. Quando o valor de

s é

menor do que

max

s , uma interação real é simulada. Após a interação (real ou delta), sorteia-

se um novo livre caminho

s e desloca-se a partícula de uma distância ),min('

max

sss = , e

assim sucessivamente. A característica Markoviana do processo de transporte garante que a

inserção de interações delta não influencia a simulação. Veja a prova disto: a probabilidade

de um passo terminar com uma interação delta é

∫

∞

−==

max

)exp( )(

sdsspp

. ( 3.31 )

Para obter a probabilidade,

dssp )( , para ocorrência da primeira interação real dentro do

intervalo ),(

dsss + , escreve-se '

max

snss

com

[

]

max

ssn

e, portanto,

max

' ss < . Esta

probabilidade é então igual à probabilidade de se ter

n interações delta sucessivas seguida

por uma interação real em uma distância dentro de )','(

dsss

da última (n-ésima)

interação delta,

dssdsspdssp

TTTT

)exp( )'exp( )(

λλλλ

−=−=

−−

, ( 3.32 )

que é o mesmo valor correto para o caso sem inserção de interações delta

(SALVAT, 2001).

3.2.6 Transporte de elétrons

No transporte de elétrons deve-se fazer distinção entre modelos de classe I e classe II. No

modelo de classe II (FIGURA 3.5), um elétron de energia inicial

E percorre uma

distância

t e então cria um elétron de energia

E . Imediatamente depois de criar este

elétron a energia do elétron primário é

EtLE

col

−−

, onde

col

L é o poder de frenamento

de colisão restrito a partículas secundárias com energia menor do que

col

tL é a

energia perdida continuamente e depositada ao longo do trajeto

t. A criação do elétron

secundário, assim como o espalhamento múltiplo, causa a mudança da direção do elétron

primário. No modelo de classe I, o elétron percorre uma distância

t e cria um elétron

secundário em algum ponto ao longo deste trajeto. A energia no final deste passo não é

explicitamente afetada pela criação deste elétron, mas é reduzida de uma quantidade

determinada de acordo com uma distribuição para perda de energia. Se a dispersão na

perda de energia fosse desconsiderada, a energia perdida seria

col

tS , o comprimento do

trajeto multiplicado pelo poder de frenamento irrestrito (excluindo efeitos radiativos). Nos

algoritmos classe I a criação de elétrons secundários não afeta explicitamente a direção do

elétron primário, contudo a deflexão associada ao espalhamento múltiplo pode levar em

conta estas deflexões.

FIGURA 3.5: Distinção entre o mecanismo de perda de energia não-correlacionado

usado em algoritmos de classe I (topo) e o mecanismo de perda de energia

correlacionado usado nos algoritmos de classe II (base).

TRANSPORTE DE ELÉTRONS

DEFINE ESTADO INICIAL DO

ELÉTRON PRIMÁRIO

RETIRA PARTÍCULA SECUNDÁRIA DO

TOPO DA PILHA COM AS

RESPECTIVAS VARIÁVEIS DE

ESTADO

E > E

corte

E ELÉTRON ESTÁ DENTRO

DO SISTEMA?

CLASSE II?

DETERMINA DISTÂNICA ATÉ A

PRÓXIMA INTERAÇÃO DISCRETA

SELECIONA O PASSO ENTRE

ESPALHAMENTOS MÚLTIPLOS

DETERMINA O ÂNGULO DE

DEFLEXÃO E MUDA A DIREÇÃO

perdida

= t L

E = E - E

perdida

O ELÉTRON SAIU DO SISTEMA?

LCANÇOU O PONTO DA INTERAÇÃO

DISCRETA?

DETERMINA O QUE OCORRE:

PRODUÇÃO DE

BREMSSTRAHLUNG OU

ELÉTRON SECUNDÁRIO

DETERMINA A ENERGIA E A

DIREÇÃO DA PARTÍCULA

SECUNDÁRIA, ARMAZENANDO

OS PARÂMETROS NA PILHA

CLASSE II?

MUDA A ENERGIA E A DIREÇÃO DA

PARTÍCULA PRIMÁRIA COMO

RESULTADO DA INTERAÇÃO

PILHA VAZIA?

TERMINA A

HISTÓRIA

DETERMINA E

perdida

FAZENDO

E = E - E

perdida

FOI CRIADA PARTÍCULA

SECUNDÁRIA?

O ELÉTRON SAIU DO SISTEMA?

SELECIONA O PASSO ENTRE

ESPALHAMENTOS MÚLTIPLOS E

PERCORRE ESTA DISTÂNCIA

DETERMINA O ÂNGULO DE

DEFLEXÃO E MUDA A DIREÇÃO

E > E

corte

E > E

corte

FIGURA 3.6: Fluxograma representando algoritmos de Monte Carlo para simulação do transporte de elétrons.

A FIGURA 3.6 mostra o fluxograma simplificado esquematizando a simulação do

transporte de elétrons. Este fluxograma mostra apenas a parte do código relativa ao

transporte; para se obter as grandezas de interesse necessita-se de um algoritmo adicional.

Outro conceito introduzido pelo fluxograma é o de “pilha” de partículas, que é

simplesmente um procedimento para se armazenar os parâmetros do espaço de fase das

partículas geradas durante a simulação de forma que estes possam ser acessados

posteriormente. Este procedimento é possível, pois este é um processo Markoviano, ou

seja, em qualquer ponto da simulação, o futuro de um fóton ou elétron independe de sua

história anterior (ROGERS, 1990).

3.2.7 Transporte de fótons

A FIGURA 3.7 apresenta um fluxograma que esquematiza a simulação do transporte

Monte Carlo de fótons. A palavra DETERMINA do diagrama simboliza um importante

conceito. Neste ponto, determina-se os parâmetros do evento a partir de distribuições de

probabilidade apropriadas.

A história de um fóton termina quando ele é absorvido ou quando sua energia se torna

menor do que a energia de corte ou quando ele deixa o volume de interesse. Os detalhes do

que fazer quando a história termina e como definir a energia de corte dependem das

informações que se deseja obter na simulação.

A FIGURA 3.7, assim como a FIGURA 3.6 para o transporte de elétrons, mostra apenas

um algoritmo para simulação do transporte do fóton, o que de fato tem pouca utilidade uma

vez que nada é registrado. Para se obter uma grandeza de interesse deve-se monitorar

certos aspectos do processo de transporte e registrar os itens relevantes (ROGERS, 1990).

3.3 Médias estatísticas e incertezas

Uma das maiores vantagens das técnicas de Monte Carlo é a possibilidade do acesso a

muito mais grandezas do que aquelas que são fisicamente mensuráveis, exceto sob

algumas raras circunstâncias. Por exemplo, utilizando-se técnicas de Monte Carlo, é

possível não só determinar a energia depositada em um tanque de água como também

conhecer a origem dos elétrons ou fótons que depositaram esta energia e quantas vezes esta

partícula ou suas descendentes sofreram um espalhamento Compton.

DEFINE ESTADO INICIAL DO FÓTON

PRIMÁRIO

TRANSPORTE DE FÓTONS

RETIRA PARTÍCULA DO TOPO DA PILHA

COM AS RESPECTIVAS VARIÁVEIS DE

ESTADO

A ENERGIA DO FÓTON É MENOR DO

QUE E

corte

DETERMINA A DISTÂNCIA ATÉ A PRÓXIMA

INTERAÇÃO

TRANSPORTA O FÓTON LEVANDO EM

CONTA A GEOMETRIA DO MEIO

O FÓTON SAIU DA REGIÃO DE

INTERESSE?

DETERMINA O QUE OCORRE:

- EFEITO FOTOELÉTRICO

- ESPALHAMENTO COMPTON

- PRODUÇÃO DE PARES

- ESPALHAMENTO COERENTE

DETERMINA A ENERGIA E A DIREÇÃO DAS

PARTÍCULAS SECUNDÁRIAS PRODUZIDAS,

ARMAZENANDO OS PARÂMETROS DESTAS

NA PILHA

PILHA VAZIA?

TERMINA A

HISTÓRIA

FIGURA 3.7: Fluxograma representando algoritmos de Monte Carlo para simulação do transporte de fótons.

Para ir direto ao ponto, considere a simulação de um feixe de elétrons de alta energia

incidindo na superfície de um fantoma de água semi-infinito. Cada elétron primário origina

um chuveiro de elétrons e fótons que são simulados um-a-um até chegarem à energia de

absorção correspondente. Para se avaliar o valor de qualquer grandeza de interesse

deve-se calcular a média desta em um grande número de chuveiros

N . Formalmente, Q

pode ser escrito como uma integral da forma (1.44)

∫

= dq qp qQ )( , ( 3.33 )

onde a FDP )(

qp é, em geral, desconhecida. A simulação de chuveiros individuais

proporciona um método prático de se sortear valores para

q a partir da FDP ‘original’:

para cada chuveiro gerado obtém-se um valor aleatório

q que é distribuído de acordo com

)(

qp . A única diferença em relação à integração de Monte Carlo considerada acima é que

neste caso a FDP )(

qp descreve uma cascata de eventos de interação, que possuem, cada

um, uma FDP característica. A estimativa de

Q é então

∑

. ( 3.34 )

Por exemplo, a energia média

dep

E depositada no fantoma por elétron incidente é dada por

∑

idep

, ( 3.35 )

onde

e é a energia depositada por todas as partículas do i-ésimo chuveiro. A incerteza

estatística (desvio padrão) para a estimativa de Monte Carlo [equação. (1.52)] é

⎥

⎦

⎤

⎢

⎣

⎡

−==

∑

NNN

)var(

. ( 3.36 )

O resultado de uma simulação deve usualmente ser escrita na forma

3± , de maneira

que o valor real de Q estará contido no intervalo )3,3(

σσ

+− com 99,7% de

probabilidade. Note que para calcular o desvio padrão (3.36) deve-se registrar os valores

q durante a simulação. Em certos casos, em que

q pode assumir apenas os valores 0 e

1, o desvio padrão pode ser determinado sem registro do quadrado destes valores:

)1(

−=

. ( 3.37 )

A simulação e o registro de grandezas durante a simulação podem também ser utilizados

para o cálculo de distribuições contínuas. O método mais simples de se fazer isto é

discretizando as distribuições, tratando-as como histogramas, e determinando as alturas de

cada uma das barras. Isto é, considere a distribuição de dose em profundidade )(

zD ,

definida como a energia média depositada por unidade de espessura do fantoma de água

por elétron incidente neste.

dzzD )( é a energia média depositada na porção do fantoma

entre

z e dzz + por elétron incidente e a integral de )(zD de 0 a ∞ é igual à energia

média depositada

dep

E (novamente, por elétron incidente). Como parte da energia

incidente é refletida pelo fantoma (através da radiação retroespalhada),

dep

E é menor do

que a energia cinética

cin

E dos elétrons incidentes. Para determinar )(zD em uma

profundidade do fantoma entre

z e

max

, deve-se, antes de qualquer coisa,

estabelecer a partição do intervalo ),0(

max

z em

caixas (bins) ),(

1 kk

−

, com

max10

0 zzzz

<<<<= K . Seja

kij

a quantidade de energia depositada na k-ésima caixa

pela

j-ésima partícula do i-ésimo chuveiro, a energia média depositada na k-ésima caixa

(por elétron incidente) é igual a

∑

kik

com

∑

≡

kijki

. ( 3.38 )

Este valor é afetado pela incerteza estatística

⎥

⎦

⎤

⎢

⎣

⎡

−=

∑

kkiEk

. ( 3.39 )

A distribuição de dose em profundidade de Monte Carlo )(

é uma função contínua

por partes,

DkkMC

DZD

3)(

para

zzz

−1

( 3.40 )

com

−

≡ e

σσ

−

≡ . ( 3.41 )

Note que a energia média em cada caixa e o desvio padrão associado têm que ser divididos

pela largura da caixa para obter a distribuição Monte Carlo final. Definida desta maneira,

)(

dá uma estimativa da dose média em cada caixa. A limitação deste método é que

aproxima-se uma distribuição contínua )(

zD por um histograma com barras de largura

finita. A princípio, seria possível obter uma melhor aproximação diminuindo-se a largura

das caixas, porém, é necessário tomar cuidado na escolha destas larguras visto que as

incertezas estatísticas podem ocultar completamente as informações no caso de caixas mais

estreitas.

Para calcular a energia média depositada e o desvio padrão associado em cada caixa,

equações (3.38) e (3.39), deve-se registrar, ao longo da simulação, os valores

,ki

e .

Porém, há casos em que o uso indiscriminado desta receita faz aumentar

consideravelmente o tempo de simulação. Considere, por exemplo, a simulação da

distribuição de dose 3D no fantoma, que pode envolver vários milhares de caixas. Para

cada caixa, as energias

kij

depositadas por cada uma das partículas do chuveiro devem ser

acumuladas em um contador parcial para obter a contribuição do chuveiro

, e, após o

término do chuveiro, o valor de

e o seu quadrado devem ser adicionados nos

contadores cumulativos. Mas como apenas uma pequena fração das caixas recebem energia

quando

um chuveiro é simulado, não é adequado tratar todas caixas da mesma maneira. O

método mais rápido para se fazer isto é transferir os registros parciais para os contadores

cumulativos apenas quando o contador parcial for receber uma contribuição de um novo

chuveiro. Isto pode ser facilmente implementado se, para cada grandeza de interesse

forem definidos três contadores reais,

Q, Q2e QP, e um indicador inteiro LQ; todos

definidos inicialmente iguais a zero. Os registros parciais

q das partículas de um chuveiro

são acumulados no contador parcial

QP, enquanto a contribuição global do chuveiro

q e o

seu quadrado são acumulados em

Q e Q2, respectivamente. A cada chuveiro é associado

um indicador – seu número,

i, dentro da seqüência de chuveiros, por exemplo – que é

registrado em

LQ quando ocorre a primeira contribuição do chuveiro para o contador

parcial

QP. No decorrer da simulação, o valor de QP é transferido para os contadores

globais

Q e Q2 apenas quando é necessário armazenar uma contribuição

q de um novo

chuveiro. Assim, o código FORTRAN para o registro de

Q ficaria:

IF(i.NE.LQ) THEN

Q=Q+QP

Q2=Q2+QP**2

QP=q

LQ=i

ELSE

QP=QP+q

ENDIF

No final da simulação, o conteúdo residual de

QP deve ser transferido para os contadores

globais.

Para algumas grandezas (como por exemplo, o número médio de eventos de espalhamento

por trajetória, a função de dose em profundidade, etc.) quase todas histórias do chuveiro

contribuem com a contagem nas caixas e as incertezas estatísticas dos resultados da

simulação se tornam relativamente pequenas. Outras grandezas (como as distribuições

angulares e de energia para partículas transmitidas através de uma lâmina espessa)

apresentam incertezas estatísticas consideráveis (isto é, grandes variâncias) porque apenas

uma pequena fração das histórias simuladas contribui para os contadores parciais

(SALVAT, 2001).

3.4 Exatidão do código

Em um código de Monte Carlo, definitivamente não importa o quanto o modelo usado é

fidedigno do ponto de vista físico ou o quão bem ele funciona, a qualidade dos resultados é

limitada pela qualidade das seções de choque utilizadas. Além dos erros intrínsecos dos

dados de seção de choque, erros adicionais têm origem no modo como o código de Monte

Carlo tabela, interpola e acessa os dados.

Ao lado dos erros nos dados de seção de choque, várias outras fontes de erro são inerentes

aos cálculos de Monte Carlo, de natureza tanto sistemática quando estatística.

3.4.1 Dados de seção de choque

A incerteza produzida no resultado final pelos erros nas seções de choque são muito

difíceis de se determinar. Veja este exemplo: um erro de 50% no poder de frenamento de

um elétron, para energias menores do que 2 MeV, poderia não ter nenhum efeito numa

curva de dose em profundidade para fótons de 10 MeV porque estes elétrons têm um

pequeno alcance e depositam toda sua energia localmente (se comparado com o poder de

penetração do fóton primário). Portanto, deve-se sempre analisar o problema em questão

para determinar como as incertezas nas seções de choque afetam a incerteza global.

Fótons

A estimativa das incertezas para as seções de choque totais para fótons vêm de Storm e

Israel (1970) e Hubbell (1982) (citados por ROGERS, 1990). Na região em que o efeito

fotoelétrico predomina (FIGURA 2.1) a incerteza é de 3 – 5%. Embora exista uma

incerteza de 3% associada ao espalhamento coerente nesta região, a incerteza global é

completamente dominada pela seção de choque fotoelétrica. Em geral, na faixa de energia

em que o efeito Compton predomina as seções de choque são precisas dentro de

aproximadamente 2 – 3%. Contudo, para elementos leves normalmente associados com

dosimetria e terapia, a incerteza pode ser menor do que 1%. A seção de choque para

produção de par deveria ser exata com uma margem de aproximadamente 5%, mas na

região de energia acima de 30 MeV deve-se analisar as incertezas relacionadas à omissão

dos processos fotonucleares.

Elétrons

Como discutido no ICRU 37 (1984) (citado por ROGERS, 1990), a incerteza no poder de

frenamento de colisão para elétrons e pósitrons é governada pela incerteza na seção de

choque do espalhamento inelástico e pelas aproximações das teorias estatísticas para a

soma das várias pequenas perdas de energia. Em baixas energias (10 – 100 keV), o ICRU

estima que a incerteza é de 2 – 3% para materiais de baixo Z e de 5 – 10% para alto Z.

Acima de 100 keV, a incerteza estimada é de 1 – 2%. A incerteza no poder de frenamento

radiativo tem origem na seção de choque de bremsstrahlung e, de acordo com o ICRU, ela

é de aproximadamente 5% para energias menores do que 2 MeV e 2 – 5% para energias

entre 2 e 50 MeV (ROGERS, 1990).

3.4.2 Erros Sistemáticos

Erros de programação

Códigos de Monte Carlo são sistemas integrados de programas computacionais grandes e

complicados. Nenhum grande código é livre de erros ou trabalha exatamente como o

planejado. A maior vantagem dos códigos que já são amplamente utilizados é que eles já

foram testados em uma variedade de aplicações.

Incertezas do Modelo

Incertezas são muitas vezes associadas com os componentes analíticos da simulação e

refletem as aproximações utilizadas. Assim como os erros inerentes do modelo empregado,

existem limitações bem conhecidas para muitos dos modelos usados. Por exemplo, muitos

códigos não incluem espalhamento coerente ou os efeitos da ligação do elétron alvo no

espalhamento Compton. Isto pode não ter nenhum efeito em muitos casos, mas poderia ser

crucial em certas aplicações. Portanto, as limitações de um código devem ser

cuidadosamente analisadas antes de utilizá-lo.

Erros de truncamento

Às vezes podem surgir dificuldades por causa da precisão finita dos números reais

representados por computadores digitais. Por exemplo, se muito mal codificado, o ponto

flutuante

de 32 bit e precisão única

∑

000.000.1

xx , ( 3.42 )

onde

−

x para todo i , dá a resposta 00.1

x e não 1.01 como esperado. Variáveis

cumulativas usadas desta forma devem ser definidas como sendo de precisão dupla. Erros

de truncamento podem causar problemas no transporte de partículas em meio com

geometria definida: uma partícula deveria estar em uma certa região e o truncamento faz

com que ela se encontre em outra. Códigos de geometria bem projetados minimizam este

tipo de erro e podem detecta-los e corrigi-los caso ocorram. Códigos de precisão dupla

ajudam neste caso, mas não eliminam o problema (ROGERS, 1990).

3.4.3 Incertezas estatísticas

Um resultado de Monte Carlo, como qualquer resultado experimental, não tem nenhum

significado a não ser que alguma estimativa de sua incerteza estatística seja dada. Como

não há uma teoria rigorosa de como determinar esta incerteza, é comum dividir uma série

de simulações Monte Carlo em um certo número

N de grupos independentes, cada um

considerando o mesmo número de histórias, e calcular um valor estimado

x para cada

grandeza de interesse em cada grupo. Com

N estimativas para

x , o valor final será igual

à média

x destes. Assumindo que

x é obtido a partir de uma distribuição normal, a

melhor estimativa da variância

é dada por [equação (3.36)]:

⎥

⎦

⎤

⎢

⎣

⎡

−=

∑

)(

. ( 3.43 )

Para

10~≥

N , pode-se assegurar que o intervalo

(

)

−

xx , contém o valor real de

com 68% de probabilidade, ou com 95% de probabilidade para s2± . Como existe uma

dependência de

com a escolha de

N (que deve ser maior ou igual a 10), deve-se

sempre indicar o valor de

N utilizado.

Método de representação onde sinaliza-se um número através de uma continuidade de algarismos e o

expoente no quadeve-se elevar o número.

Para combinar os resultados e incertezas de

m simulações independentes com valores

distintos para

N , utiliza-se as seguintes fórmulas:

∑

( 3.44 )

∑

⎟

⎠

⎞

⎜

⎝

⎛

σσ

( 3.45 )

onde

N é o número total de histórias e o índice k indica que aquela variável é associada à

k-ésima simulação. Incertezas maiores ou iguais a 5% tendem a ser ligeiramente

subestimadas e incertezas menores do que isto parecem ser mais precisas e concordam

qualitativamente com o que seria esperado se a distribuição de

x fosse realmente normal

(ROGERS, 1990).

3.4.4 Eficiência Computacional

Assumindo que as estimativas das várias simulações independentes são extraídas de uma

distribuição normal, no limite de um grande número de histórias a grandeza

N é uma

constante. Assim, em um determinado cálculo, a incerteza estimada é aproximadamente

proporcional ao inverso da raiz quadrada do número de histórias utilizado. Esta relação

pode ser usada para determinar o número de histórias necessário para se alcançar um

resultado com uma certa precisão, desde que já exista uma estimativa, mesmo que

estatisticamente pobre, para esta incerteza.

A relação mencionada acima permite também uma definição prática para a eficiência de

um cálculo de Monte Carlo:

σε

T= ( 3.46 )

onde

é o tempo computacional necessário para obter uma variância igual a

. Se a

técnica para os cálculos é mantida a mesma,

permanece constante mesmo aumentando-

se o número de histórias

N . A definição desta grandeza é útil porque quando se

introduzem técnicas de redução de variância, o tempo de computação para um determinado

número de histórias irá, em geral, aumentar. Porém, se o produto

T diminui, há uma

melhoria na eficiência computacional (o tempo necessário para alcançar um dado valor

torna-se menor ou, alternativamente, para um dado tempo, o valor obtido para

reduzido) (ROGERS, 1990).

3.5 Redução de Variância

No transporte Monte Carlo de elétrons e fótons gasta-se muito tempo gerando números

aleatórios, girando os sistemas de coordenadas, determinando o ângulo de espalhamento

múltiplo, traçando a geometria, etc. Um código de Monte Carlo com uma codificação

eficiente e uma boa estrutura pode diminuir consideravelmente o tempo computacional em

uma simulação. Por exemplo, pode-se tornar os cálculos mais rápidos construindo-se

tabelas para funções mais freqüentemente usadas, assim como melhorar a eficiência dos

procedimentos de amostragem de processos físicos. Por outro lado, técnicas de redução de

variância podem melhorar a eficiência dando mais ênfase às grandezas físicas de interesse,

produzindo então uma maior quantidade de informação relevante em um determinado

tempo computacional (ROGERS, 1990).

A princípio, o erro estatístico de uma grandeza pode ser de certa forma reduzido (sem

aumentar o tempo computacional de simulação) através do uso de técnicas de redução de

variância. Infelizmente, estas técnicas de otimização são extremamente dependentes do

problema a ser tratado, o que impossibilita a definição de um procedimento geral para

minimizar a variância. Por outro lado, não se deve superestimar a importância da redução

de variância, pois em muitos casos, a simulação análoga

pode ser realizada em um tempo

razoável. Gastar tempo complicando o programa para se obter apenas uma modesta

redução no tempo computacional pode não ser um bom investimento. É também

importante compreender que um método de redução de variância eficiente usualmente

diminui o erro estatístico de uma dada grandeza

Q às custas do aumento das incertezas de

outras grandezas. Portanto, as técnicas de redução de variância não são recomendadas

quando se busca uma descrição global do processo de transporte (SALVAT, 2001). Nesta

seção é feita uma breve descrição de algumas das técnicas que, com uma programação

simples, podem ser úteis para melhorar a solução de alguns problemas.

3.5.1 Interação Forçada

Uma grande variância pode ser resultado de probabilidades de interação extremamente

baixas. Um exemplo disto é a simulação do espectro de energia de fótons de

bremsstrahlung emitidos por elétrons de energia intermediária (~ 100 keV) em uma lâmina

fina de um determinado material. Como eventos radiativos são muito menos prováveis do

Aqui, o termo “análoga” refere-se a simulações detalhadas, condensadas ou mistas, que não incorporam

procedimentos de redução de variância.

que os espalhamentos elástico e inelástico, a incerteza do espectro do fóton simulado será

relativamente grande. Nestes casos, um eficiente método para redução da variância é o

aumento artificial da probabilidade de interação para o evento “A” de interesse. A

implementação prática da interação forçada consiste da substituição do livre caminho

médio

do processo real por um menor,

fA,

, forçando as interações do tipo A

ocorrerem mais freqüentemente do que no processo real. Considera-se que as FDP’s para

perda de energia e deflexões angulares (e as direções de emissão das partículas

secundárias, se for o caso) nas interações forçadas são as mesmas das interações reais e

para se sortear a distância a ser percorrida até a próxima interação faz-se uso da

distribuição exponencial com o livre caminho médio reduzido

fA,

. Isto é equivalente a

aumentar a probabilidade de interação, através do processo A, por unidade de distância

percorrida, pelo fator

. ( 3.47 )

Para que esta medida não interferira nos resultados da simulação, deve-se corrigir as

distorções introduzidas como a seguir:

(i) Um peso

)1(

w é associado a cada partícula primária. As partículas secundárias

produzidas em interações forçadas têm peso F

)1()2(

ww = , e os pesos das

sucessivas gerações de partículas secundárias “forçadas” são F

)1()( −

ww . Às

partículas secundárias geradas em interações não forçadas (geradas por outros

processos que não o A) são dados pesos iguais ao da sua partícula mãe.

(ii) O peso

)()( k

ww = é dado à deposição de energia (e qualquer outra alteração

do meio, como por exemplo, deposição de carga) resultante de interações forçadas

de uma partícula de peso

)(k

w . Para interações não forçadas,

)()( k

ww = .

(iii) Quando interações forçadas são simuladas para se determinar a perda de energia e a

possível emissão de radiação secundária, as variáveis de estado da partícula em

questão é alterada com a probabilidade

F1 . Ou seja, a energia

e a direção de

movimento

do projétil somente variam quando o valor

do número aleatório

for menor do que

F1 , caso contrário,

e d

permanecem inalterados.

Considere um problema de blindagem em que, na média, apenas um fóton dentre mil

atravessa o bloqueio. Mudando a seção de choque de maneira que, para cada 1000

histórias, 10 fótons atravessem a blindagem, a eficiência aumentará, uma vez que a

incerteza percentual no número de fótons atravessando o bloqueio será reduzida de

) 11 ( %100 ±± para ) 1010 ( %32 ±± . Neste caso, deve-se ter em mente que a

probabilidade de transmissão foi multiplicada por 10 e, assim, a resposta deve ser dividida

pelo mesmo fator, ficando 3,00,1 ±

w e

q simbolizam o peso e a contribuição da i-ésima partícula primária no cálculo de

alguma grandeza

Q, e

w e

q (1>j ) representam os pesos e contribuições das j-ésimas

partículas secundárias geradas pela

i-ésima primária, a estimativa Monte Carlo de Q

obtida a partir da simulação de

N histórias será

∑

ijij

. ( 3.48 )

Evidentemente, as estimativas de Q obtidas com interações forçadas e com uma

simulação análoga são iguais (do ponto de vista estatístico, ou seja, no limite

∞

→N , a

diferença entre elas tende a zero). O desvio padrão será dado por

⎥

⎦

⎤

⎢

⎣

⎡

−

⎟

⎠

⎞

⎜

⎝

⎛

∑∑

ijijQ

Qqw

. ( 3.49 )

As grandezas diretamente relacionadas às interações forçadas terão erro estatístico

reduzido devido ao aumento no número destas interações, porém, para um dado tempo de

simulação, outras grandezas podem apresentar desvios padrões maiores do que aqueles da

simulação análoga em função do tempo gasto simulando-se interações forçadas

(SALVAT, 2001).

3.5.2 Cisão (splitting) e roleta Russa

Estas duas técnicas, que são normalmente usadas em conjunto, são eficientes em

problemas cujo interesse está focado em uma determinada região do espaço. O cálculo de

funções de dose em profundidade em um objeto irradiado e, no caso de feixes de radiação

colimados, o cálculo de isodoses em pontos distantes do eixo do feixe são exemplos típicos

de problemas tratáveis com estas técnicas. A idéia básica dos métodos de cisão e roleta

Russa é favorecer o fluxo de radiação para a região de interesse e impedir que a radiação

deixe esta região. Estas técnicas são úteis também em problemas em que é necessária

apenas uma descrição parcial do processo de transporte. A “região de interesse” pode então

ser um volume limitado no espaço das variáveis de estado (

,,E ); assim, em estudos de

retroespalhamento de radiação, a região de interesse pode ser estabelecida como sendo a

região da amostra mais próxima à superfície irradiada e o conjunto de partículas cuja

direção aponta para esta superfície.

Como no caso da interação forçada, aqui, a redução de variância é realizada modificando-

se os pesos das partículas. Assume-se que as partículas primárias começam com peso

unitário e a cada partícula secundária produzida por uma primária é atribuído um peso

inicial igual ao desta primária. No método de cisão, uma partícula, com peso

w e em um

determinado estado, é transformada em um número

de partículas idênticas com

pesos

ww = no mesmo estado inicial. A cisão deve ser aplicada quando a partícula se

aproxima da região de interesse. A técnica da roleta Russa é, de certo modo, o contrário da

cisão: quando a partícula tende a afastar-se da região de interesse, ela pode ser “destruída”

com uma probabilidade 1<

, e, se isto não ocorrer, seu peso é aumentado de um fator

)1(1

− . Aqui, destruir a partícula significa simplesmente descarta-la, não considerando

mais suas contribuições para os cálculos. Estes dois métodos não influenciam os resultados

da simulação, sendo a média e o desvio padrão das grandezas calculadas dadas pelas

equações (3.48) e (3.49). A eficiência destes métodos está relacionada com os valores

adotados para os parâmetros

, e no critério usado para decidir quando a cisão e a

eliminação da partícula devem ser aplicadas (SALVAT, 2001).

3.5.3 Rejeição de alcance

O método de rejeição de alcance consiste em absorver uma partícula quando não houver

possibilidade desta partícula (e suas possíveis secundárias) alcançar ou deixar a região de

interesse. Este método pode ser útil, por exemplo, quando se deseja calcular a energia total

depositada por elétrons (ou pósitrons) em uma certa região do espaço. Quando o alcance

residual de uma partícula for menor do que a distância até a mais próxima superfície que

delimita a região de interesse, considera-se que toda sua energia é depositada dentro ou

fora desta região (dependendo se ela estava no interior desta ou não) e a simulação da

trajetória pode ser interrompida (SALVAT, 2001). Em contrapartida, este método introduz

uma aproximação ao ignorar a possibilidade do elétron produzir um fóton que poderia

escapar ou entrar na região de interesse; por esta razão usualmente determina-se uma

energia máxima para aplicação deste método (ROGERS, 1990). Além disto, a rejeição de

alcance não é adequada para o transporte de fótons, uma vez que o conceito de alcance

para fótons não é bem definido (ou, mais precisamente, as flutuações no tamanho da

trajetória para fótons são muito grandes) (SALVAT, 2001).

3.5.4 Outros métodos

Freqüentemente, uma redução de variância efetiva pode ser obtida simplesmente evitando-

se cálculos desnecessários; como geralmente acontece nos códigos de simulação que

incorporam pacotes de geometria genéricos. Para geometrias simples (por exemplo, planar,

esférica ou cilíndrica) o programa pode ser substancialmente simplificado e isto pode

acelerar a simulação consideravelmente. De maneira geral, o uso adequado de possíveis

simetrias do problema em consideração pode reduzir a variância apreciavelmente

(SALVAT, 2001).

METODOLOGIA

O problema tratado neste trabalho consiste na simulação do transporte de elétrons de baixa

energia (110 eV a 50 keV), e dos elétrons e fótons secundários gerados por estes,

provenientes de uma fonte monoenergética isotrópica localizada no interior de um meio

tecido equivalente homogêneo semi-infinito (suas dimensões são muito maiores do que o

alcance das partículas em questão). A composição deste meio é igual à definida pelo ICRP

(International Commission on Radiation Protection) como tecido mole (ver apêndice A).

As simulações foram realizadas utilizando-se o código de Monte Carlo PENELOPE que

será descrito a seguir. A escolha deste código em especial para realização do estudo do

alcance de elétrons de baixa energia se deu devido ao fato dele ser capaz de simular o

transporte de elétrons até energias tão baixas quanto 100 eV, ao passo que códigos como o

EGS4, ETRAN, ITS3 e MCNP são, de maneira geral, limitados a elétrons com energias

acima de 1 a 20 keV.

Foram realizadas, separadamente para cada energia, as simulações do transporte de

elétrons com energias iniciais de 110, 200, 300, 400, 500, 600, 700, 800 e 900 eV, e 1, 5,

10, 25 e 50 keV (o valor de 110 eV foi escolhido devido à limitação do código de Monte

Carlo utilizado que simula elétrons somente até o limite inferior de 100 eV). Cada

simulação era composta por um certo número N de chuveiros (partículas primárias)

definido de acordo com a energia inicial dos elétrons primários emitidos pela fonte. Esta

diferença no número de chuveiros simulados para cada energia se dá devido ao fato de que

quanto maior a energia inicial do elétron, maior o número de interações e partículas

secundárias produzidas até a degradação completa da energia deste elétron e, portanto,

maior o número de informações a serem processadas e armazenadas. Assim, para cada

energia, o número de chuveiros é limitado em função do esforço e memória computacional

necessários. Portanto, quanto maior a energia inicial do elétron, menor o número de

chuveiros em cada simulação. A determinação de N foi feita testando-se o limite em que o

código PENELOPE e o próprio computador era capaz de trabalhar, simulando um número

cada vez maior de chuverios para cada energia. Finalmente, para cada energia considerada,

foram realizadas dez simulações, contendo cada uma N chuveiros.

Durante a simulação, as informações a respeito de toda história de cada partícula simulada,

como posição em que ocorreram interações, tipo de interação, energia depositada em cada

evento, direção de deflexão, etc., foram gravadas seqüencialmente em arquivos sem

formato (binário) para permitir o acesso a estas informações posteriormente. Os dados

foram gravados na representação binária porque o grande número de interações e partículas

secundárias geradas no transporte de elétrons produz um volume tão grande de informação

que não é possível armazená-las da maneira convencional em arquivos de texto, que

inclusive permitiriam a visualização direta destas informações. No PENELOPE, cada

partícula é simulada até que sua energia se torne menor do que a energia de absorção desta

no meio. A simulação inicia com uma partícula primária e após o término da história desta,

são simuladas as histórias das partículas secundárias geradas, uma de cada vez,

completando o chuveiro. Terminadas as partículas secundárias, é iniciado um novo

chuveiro e assim sucessivamente. Desta forma, o arquivo em que foram gravadas as

informações relevantes das simulações possui a seguinte estrutura:

(i) Na primeira linha, o estado inicial (energia, posição e direção iniciais) da partícula

primária a ser simulada;

(ii) A cada interação sofrida por esta partícula primária uma nova linha é gravada

(contendo o estado desta partícula após a interação, o mecanismo de interação que

esta sofreu, a energia perdida, etc.) repetidamente até o término da história desta

(que ocorre quando sua energia é completamente degradada);

(iii) Terminada a história da partícula primária, na linha seguinte é gravado o estado

inicial da primeira partícula secundária a ser simulada, sendo gravadas nas linhas

subsequentes as informações relativas às interações sofridas por esta até o término

de sua história. Este passo se repete até a simulação de todas as secundárias

finalizando assim um chuveiro;

(iv) Terminadas as histórias das secundárias, é iniciado um novo chuveiro e a

gravação prossegue repetindo todas as etapas a partir do passo (i) até a simulação

do número N de chuveiros (partículas primárias) pré-determinado.

A partir destes dados obtidos nas simulações (posição em que ocorreram as interações, tipo

de interação ocorrida, ângulo de deflexão, energia transferida, etc.) foram então calculados

os alcances para elétrons com as energias iniciais mencionadas anteriormente.

4 O código PENELOPE

O PENELOPE é um algoritmo de Monte Carlo para a simulação do transporte acoplado de

elétrons e fótons. Este nome é uma abreviatura de “PENetration and Energy Loss of

Positrons and Electrons” – penetração e perda de energia de elétrons e pósitrons (a

simulação de fótons foi introduzida posteriormente). O algoritmo de simulação é baseado

em um modelo de espalhamento que combina bases de dados numéricas com modelos

analíticos de seções de choque para os diferentes mecanismos de interação e trabalha com

energias (energia cinética no caso de elétrons e pósitrons) entre aproximadamente 100 eV e

1 GeV. Os processos de interação relevantes para fótons são: espalhamento coerente

(Rayleigh), espalhamento incoerente (Compton), efeito fotoelétrico e produção do par

elétron-pósitron. Outras interações (como, por exemplo, absorção fotonuclear) ocorrem

com uma probabilidade muito menor e podem ser desconsideradas na maior parte dos

casos práticos. Para elétrons e pósitrons, as interações consideradas são: espalhamento

elástico, colisões inelásticas e emissão de bremsstrahlung; pósitrons podem também sofrer

aniquilação, tanto em repouso quanto em vôo. O código PENELOPE é particularmente útil

em aplicações que envolvem o transporte de radiação em outros materiais que não água

(STEWART, 2002).

4.1 Interações de fótons

4.1.1 Espalhamento coerente (Rayleigh)

Espalhamento coerente ou Rayleigh é o processo pelo qual fótons são espalhados por

elétrons atômicos (ligados) sem excitação do átomo alvo, isto é, a energia dos fótons

incidente e espalhado é a mesma. O espalhamento é qualificado como “coerente” porque

ele surge da interferência entre ondas eletromagnéticas secundárias vindas de diferentes

partes da distribuição de cargas atômicas.

4.1.2 Efeito fotoelétrico

No efeito fotoelétrico, um fóton de energia

é absorvido pelo átomo alvo, que faz uma

transição para um estado excitado. Os feixes de fótons encontrados em estudos de

transporte de radiação têm densidade de fótons relativamente baixa e, como conseqüência,

é observada apenas a absorção de fótons individuais (em feixes intensos de fótons de baixa

energia, como os de lasers de alta potência, é possível a absorção simultânea de vários

fótons). Para representar os estados atômicos, pode-se adotar um modelo elétron-

independente, assim como o modelo de Dirac-Hartree-Fock-Slater, em que cada elétron

ocupa um orbital, com energia de ionização bem definida. O conjunto de orbitais com os

mesmos números quânticos para o momento principal e total e a mesma paridade

constituem uma camada. Cada camada i pode acomodar um número finito de elétrons, com

energia de ionização característica

U .

A interação com o campo do fóton é considerado como uma perturbação de primeira

ordem (o que é apropriado para campos com baixa densidade de fótons) sendo permitidas

apenas transições de um único elétron; isto é, no efeito fotoelétrico, o fóton é absorvido por

um elétron individual da camada “ativa” i, que deixa o átomo com energia cinética

UEE −= . Evidentemente, a fotoionização de uma dada camada somente é possível

quando a energia do fóton incidente excede a energia de ionização correspondente; isto dá

origem aos “cortes” (edges) de absorção característicos na seção de choque fotoelétrica.

4.1.3 Espalhamento incoerente (Compton)

No espalhamento Compton, um fóton de energia

interage com um elétron atômico, que

o absorve e reemite um fóton (Compton) secundário de energia '

na direção ),(

Ω

relativa à direção do fóton incidente. Após uma interação Compton com a i-ésima camada,

o elétron alvo é ejetado com energia cinética 0 ' >

−

UEEE , onde

U é a energia de

ionização da camada em questão, e o átomo residual fica então em um estado excitado com

uma vacância na camada i.

4.1.4 Produção de pares

Pares de elétrons e pósitrons podem ser criados pela absorção de um fóton na vizinhança

de uma partícula massiva, um núcleo ou um elétron, que absorve energia e momento de

forma tal que estas duas grandezas são conservadas. A energia mínima para produção de

pares no campo de um núcleo (supostamente de massa infinita) é

2 cm

. Quando a

produção de pares ocorre no campo de um elétron, este elétron alvo, após o evento, recua

com uma energia cinética considerável; este processo é conhecido como produção de

“tripleto”. Se o elétron alvo estiver em repouso, a produção de tripleto é possível apenas

para fótons com energia maior do que

4 cm

(SALVAT, 2001).

4.2 Interações de elétrons e pósitrons

4.2.1 Colisão elástica

Por definição, interações elásticas são aquelas em que os estados quânticos inicial e final

do átomo alvo são os mesmos, em geral o estado fundamental. As deflexões angulares nas

trajetórias de elétrons na matéria são principalmente (mas não unicamente) devido ao

espalhamento elástico. Note que há uma certa transferência de energia do projétil para o

alvo, que causa o recuo deste último. Devido à grande massa do alvo (

Zm3500~), a

energia média perdida pelo projétil é uma fração muito pequena de sua energia inicial e é,

em geral, desprezada, o que é equivalente a assumir que o alvo tem massa infinita e não

recua.

4.2.2 Colisão inelástica

A colisão inelástica é o mecanismo dominante através do qual elétrons e pósitrons, de

energias baixas e intermediárias, perdem energia para o meio, produzindo excitações

eletrônicas e ionizações neste. A teoria quântica do espalhamento inelástico de partículas

carregadas por átomos individuais e moléculas foi inicialmente formulada por Bethe (1930,

1932) (citado por SALVAT, 2001) com base na aproximação de primeira ordem (onda

plana) de Born. A extensão da teoria para colisões inelásticas na matéria condensada foi

discutida por Fano (1963) (citado por SALVAT, 2001). Os aspectos formais desta teoria

são bastante complicados mas, felizmente, os resultados são essencialmente equivalentes

àqueles da teoria dielétrica clássica.

O efeito de colisões inelásticas individuais no projétil é completamente especificado pela

perda de energia

W e pelos ângulos de espalhamento polar e azimutal

respectivamente. Para meios amorfos, com átomos (ou moléculas) aleatoriamente

orientados, a SCD para colisões inelásticas é independente do ângulo de espalhamento

azimutal

4.2.3 Emissão de bremsstrahlung

Como resultado da aceleração causada pelo campo eletrostático de átomos, elétrons (ou

pósitrons) rápidos emitem fótons de bremsstrahlung (radiação de frenamento). Em cada

evento de emissão de bremsstrahlung, um elétron com energia cinética

gera um fóton de

energia

W , que assume valores no intervalo entre 0 e

. O processo é descrito por uma

SCD atômica, diferencial em relação à perda de energia

W e às direções finais do projétil

e do fóton emitido. As deflexões angulares do projétil são obtidas da SCD para

espalhamento elástico e, consequentemente, a direção de movimento do projétil permanece

inalterada na simulação de eventos radiativos.

4.2.4 Aniquilação de pósitrons

Pósitrons penetrando um meio de número atômico Z com energia cinética

podem se

aniquilar com os elétrons do meio emitindo dois fótons. Se for assumido que os elétrons

alvo estão livres e em repouso, desconsiderando então os efeitos de ligação, é possível a

aniquilação com emissão de apenas um fóton. Quando a aniquilação acontece em vôo, ou

seja, quando a energia cinética

do pósitron é maior do que a energia de absorção, os

dois fótons podem ter energias diferentes,

−

E e

E , que somadas dão

2 cmE

(SALVAT, 2001).

4.3 Relaxação atômica

O PENELOPE simula a emissão de radiação característica e elétrons Auger que resultam

de vacâncias produzidas nas camadas K e subcamadas L por absorção fotoelétrica,

espalhamento Compton e impacto com elétron/pósitron. A relaxação destas vacâncias é

seguida até que as camadas K e L sejam preenchidas, isto é, até que as vacâncias tenham

migrado para a camada M ou outra camada mais externa. As vacâncias nestas camadas

mais externas dão origem a radiações secundárias com energias muito menores, cujo

principal efeito é o de difundir a energia de excitação do íon no material em torno deste.

Para se ter uma descrição confiável da distribuição de dose, e outras características

macroscópicas do transporte, é necessário seguir apenas a radiação secundária que é capaz

de propagar distâncias da ordem de, por exemplo, 1% da distância que a partícula primária

penetra no meio (ou do alcance).

Para simplificar a descrição dos processos de ionização de camadas mais externas, assume-

se que, quando a ionização ocorre na camada M ou outra mais externa, um elétron

secundário (radiação delta) é ejetado do íon com uma energia cinética

E igual à energia

depositada pela partícula primária,

⎪

⎩

⎪

⎨

⎧

−

sitron.elétron/pó com impacto no

ca,fotoelétri absorção na

,Compton toespalhamen no '

dep

( 4.1 )

Isto é, o elétron fica com toda energia de excitação do íon e nenhuma radiação fluorescente

é simulada. Na realidade, os elétrons emitidos têm energia menor do que os valores (4.1) e

podem ser seguidos de raios-X característicos, que têm livre caminho médio geralmente

muito maior do que o alcance de Bethe para fotoelétrons. Aumentando artificialmente a

energia inicial do elétron permite-se que ele transporte energia para mais longe do íon

compensando as outras radiações desconsideradas durante a cascata de desexcitação.

No caso da ionização de uma camada interna i, isto é uma camada K ou L, considera-se

que o elétron é emitido com energia cinética

ideps

UEE

−

, ( 4.2 )

onde

U é a energia de ionização da camada ativa, ficando o átomo alvo com uma

vacância na camada i. Como mencionado, considera-se apenas raios-X característicos e

elétrons Auger emitidos nos primeiros estágios do processo de relaxação. Além disto,

assume-se que estas radiações secundárias são emitidas isotropicamente pelo átomo

excitado (SALVAT, 2001).

4.4 Mecânica do transporte de elétrons e pósitrons

No PENELOPE, o transporte de fótons é simulado por meio do método detalhado

convencional, isto é, todos os eventos de interação na história do fóton são simulados

cronologicamente. A simulação do transporte de elétrons e pósitrons é feita por meio de

um modelo classe II: interações catastróficas, com ângulo de espalhamento

e/ou perdas

de energia

W maiores do que os valores de corte

W , são simuladas detalhadamente

enquanto, interações suaves, com ângulo de espalhamento ou perda de energia menor do

que os valores de corte correspondentes, são descritas por meio de aproximações de

espalhamento múltiplo. Este método de simulação lida de forma apropriada com

deslocamentos laterais e transposição de interfaces e dá uma descrição consistente da

dispersão de energia (“energy straggling”). A simulação é estável sob variações dos

valores de corte

W , que podem assumir valores consideravelmente grandes,

diminuindo o tempo de simulação sem alterar os resultados.

Partículas secundárias emitidas com energia inicial maior do que a energia de absorção

abs

E são armazenadas e simuladas após o término da história da partícula primária. Neste

trabalho, E

abs

= 100 eV. As partículas secundárias são produzidas diretamente (em

interações catastróficas como colisões inelásticas, emissão de bremsstrahlung, aniquilação

de pósitron, espalhamento Compton, absorção fotoelétrica e produção de pares) ou

indiretamente como radiação fluorescente (raios-X característicos e elétrons Auger). O

PENELOPE simula a emissão de raios-X característicos e elétrons Auger que resultam de

vacâncias produzidas nas camadas K e subcamadas L, sendo a relaxação destas vacâncias

seguida até o preenchimento das camadas K e L, ou seja, até as vacâncias migrarem para a

camada M ou outra camada mais externa (SALVAT, 2001).

4.5 Simulação mista (classe II) do espalhamento elástico

4.5.1 Espalhamento múltiplo. Método da dobradiça aleatória (“random hinge”)

No PENELOPE, a deflexão angular e o deslocamento lateral devido as múltiplas colisões

suaves em um passo de comprimento

s são simulados através do método da dobradiça

aleatória (FERNÁNDEZ-VAREA et al., 1993 citado por SALVAT, 2001). O algoritmo

deste método pode ser programado como a seguir (FIGURA 4.1):

(i) Primeiramente, o elétron percorre uma distância aleatória

, sorteada

uniformemente no intervalo ),0( s , na direção inicial.

(ii) Ocorre então um evento de espalhamento suave artificial (uma dobradiça), em que

o elétron muda sua direção de movimento de acordo com a distribuição para

espalhamento múltiplo );(

)(

, onde

é o ângulo polar que define a nova

direção de movimento.

(iii) Finalmente, o elétron percorre uma distância

−

s nesta nova direção.

Obviamente, este algoritmo fornece a distribuição angular exata no final do passo.

4.5.2 Perdas de energia suaves. “Continuous-Slowing-Down-Approximation”

O uso do CSDA (“Continuous-Slowing-Down-Approximation”) na descrição de interações

suaves é satisfatoriamente justificável quando a dispersão da energia decorrente a estas

FIGURA 4.1: Simulação do efeito global de interações suaves entre duas colisões

catastróficas consecutivas através do método da dobradiça aleatória.

interações é desprezível. Este é o caso que ocorre, por exemplo, quando as energias de

corte

W e

W são pequenas, condição que causa também a redução da fração do poder

de frenamento relativa às interações suaves. Para melhorar a descrição de dispersão da

energia seria necessário reduzir as energias de corte, mas isto faz aumentar o número de

eventos catastróficos inelásticos e radiativos a serem simulados ao longo de cada trajetória,

aumentando assim o tempo de simulação. A proposta do PENELOPE é ir além do CSDA

introduzindo-se a dispersão da energia na descrição das interações suaves. Fica claro que,

procedendo desta maneira, será possível o uso de valores maiores para as energias de corte

W e

W , aumentando assim a velocidade da simulação sem alterar as distribuições de

energia. Em cálculos de distribuição espacial de dose, a perda de energia

w decorrente de

interações suaves pode ser considerada como sendo localmente depositada em uma posição

aleatória uniformemente distribuída ao longo do passo. Este procedimento produz

distribuições idênticas àquelas obtidas assumindo-se que a energia é depositada a uma taxa

constante ao longo do passo, porém é computacionalmente mais simples. Com base nisto,

o PENELOPE simula o efeito combinado de todas as colisões elásticas suaves que ocorrem

entre dois eventos catastróficos sucessivos, separados de uma distância

s , como um único

evento (uma dobradiça) em que a partícula muda sua direção de movimento e perde

energia de acordo com as respectivas distribuições de probabilidade. A posição da

dobradiça é uniformemente distribuída ao longo do passo, como no caso do espalhamento

elástico puro.

4.5.3 Combinação de espalhamento e perda de energia

O espalhamento múltiplo e a perda de energia suave foram tratados como processos

independentes, mas na verdade eles coexistem. É necessário, então, considerar esta

reciprocidade e definir a base de um algoritmo que simula estes efeitos combinados.

Na ausência de perdas de energia suaves, a FDP para comprimento do passo

s entre dois

eventos catastróficos sucessivos (ou de um dado ponto da trajetória até o próximo evento

catastrófico) é

)exp()( ssp

−

( 4.3 )

onde

Σ é o inverso do livre caminho médio

)(h

entre dois eventos catastróficos

consecutivos

. Em cada evento catastrófico, uma e apenas uma interação – i = “el”

(elástica), “in” (inelástica), “br” (bremsstrahlung) ou “na” (aniquilação de pósitron) –

ocorre com probabilidade

)()( h

σσ

= . ( 4.4 )

Quando perdas de energia suaves são consideradas, a FDP para a distância

s percorrida

pela partícula até a interação catastrófica seguinte não é dada pela equação (4.3), porque o

livre caminho médio

)(h

varia com a energia e pode mudar consideravelmente ao longo

de um único passo. A maneira mais simples de lidar com este problema é limitar o

comprimento do passo para garantir que a perda de energia média seja muito menor do que

a energia cinética

do início deste, e considerar que )(

)(

permanece essencialmente

constante ao longo do passo. Assim, a energia média perdida em um passo será dada por

)(

ESE

=∆ , ( 4.5 )

onde

)()()( ESESES

brin

( 4.6 )

é o poder de frenamento total. Uma vez que o livre caminho médio entre eventos

catastróficos consecutivos de qualquer tipo é menor do que o livre caminho médio entre

eventos elásticos catastróficos, a perda de energia por passo pode ser limitada redefinindo-

se o livre caminho médio catastrófico. Para permitir frações de perda de energia

EE∆ , ao

longo do passo, da ordem de

C (um pequeno valor, digamos que 0,05), basta fazer

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

)(

),(min),(max)(

21,1

)(

CECEE

elel

λλλ

, ( 4.7 )

onde

1,el

é o primeiro livre caminho médio de transporte

Isto limita efetivamente a perda de energia média por passo, porém aumenta a freqüência

dos eventos catastróficos elásticos. Os parâmetros

C e

C na equação (4.7), definidos

pelo usuário, influenciam no tempo computacional necessário para simular cada trajetória.

[]

(

)

[]

)()()()()(

)( h

σσσσσλ

+++==≡Σ

−

, onde N é o número de centros espalhadores

por unidade de volume,

)(h

é a seção de choque por átomo total para interações catastróficas.

O l-ésimo livre caminho médio de transporte é definido por

1,1,

elel

≡

, onde N é o número de átomos

por unidade de volume,

()

[]

∫

Ω

−≡ d

cos1

11,

θσ

Ω

d é o elemento de ângulo sólido na

direção

()

, e

(

)

cos

P são os polinômios de Legendre.

Em condições ideais, eles não devem ter qualquer influência na precisão dos resultados da

simulação, o que ocorre somente quando seus valores são suficientemente pequenos.

As deflexões angulares devido às interações suaves ao longo de um passo de comprimento

s são geradas a partir de uma distribuição artificial que é função do livre caminho médio

)(s

comb

, correspondente à combinação dos processos de espalhamento e perda de energia na

dobradiça. Para levar em conta (pelo menos parcialmente) a dependência energética destas

grandezas, a perda de energia e a deflexão angular (que ocorrem na dobradiça) são

consideradas como sendo processos independentes e são simuladas em ordem aleatória. Ou

seja, a deflexão angular suave é calculada considerando-se a energia tanto no começo

quanto no final do passo, com probabilidades iguais. Isto é equivalente a assumir que o

livre caminho médio

)(s

comb

varia linearmente com a energia. Este método é bastante preciso

e computacionalmente simples e tem apenas a exigência da perda fracional de energia em

cada passo (que é da ordem de

C ) ser suficientemente pequena (SALVAT, 2001).

4.6 Geração das histórias

A simulação da história de um elétron ou de um fóton é basicamente uma sucessão

cronológica de eventos. Estes podem ser eventos catastróficos, interações suaves

(dobradiças) ou outro estágio relevante da história da partícula (seu estado inicial, a

passagem através de uma interface ou sua absorção no final da trajetória). O trajeto da

partícula entre dois eventos sucessivos é uma linha reta e é chamado segmento. A porção

do trajeto entre dois eventos catastróficos é chamada passo e consiste de dois segmentos

com uma dobradiça entre eles (no caso de simulação mista).

No PENELOPE, o algoritmo de transporte de elétrons e pósitrons é controlado pelos

seguintes parâmetros de simulação:

–

C : Deflexão angular média,

cos1− , produzida por espalhamento elástico

múltiplo ao longo do livre caminho médio

)(h

entre dois eventos catastróficos de

espalhamento elástico. O ideal é que

C seja da ordem de 0,05, mas ele pode

assumir qualquer valor entre 0 (simulação detalhada) e 0,2, que corresponde a uma

deflexão angular média

º37~

após um passo de comprimento

)(h

–

C : Valor máximo da perda de energia média (percentual) entre eventos

catastróficos consecutivos de espalhamento elástico. Pode assumir valores entre 0 e

0,2, mas em geral, um valor da ordem de 0,05 é o mais adequado.

–

W : Valor de corte para perda de energia para colisões inelásticas catastróficas.

–

W : Valor de corte para perda de energia em emissões de bremsstrahlung

catastróficas.

Estes parâmetros determinam a precisão e a velocidade da simulação. Para garantir boa

precisão,

C e

C devem ter valores pequenos (da ordem de 0,01 ou menos). Com valores

grandes para estes parâmetros a simulação é mais rápida, porém ocorre uma certa perda de

precisão. As energias de corte

W e

W influenciam principalmente as distribuições de

energia simuladas. A simulação é mais rápida para energias de corte maiores, mas se estas

são muito grandes, as distribuições simuladas podem se tornar um pouco distorcidas. Na

prática, as distribuições de energia simuladas costumam ser completamente insensíveis aos

valores de

W e

W adotados quando estes são menores do que a largura das ‘caixas’

utilizadas para o cálculo das distribuições de energia. Portanto, é a resolução desejada que

determina os valores máximos para as energias de corte.

O efeito combinado de todas as interações elásticas suaves e interações para interrupção da

trajetória da partícula é simulado como um único evento artificial (dobradiça), em que a

partícula muda sua direção de movimento e perde energia. Quando

W é menor do que a

mais baixa energia de ressonância, a simulação de colisões inelásticas se torna puramente

detalhada, isto é, as colisões inelásticas não contribuem para o “poder de frenamento

suave”. Por outro lado, a simulação de emissão de bremsstrahlung é possível apenas

através do esquema misto devido à divergência da SCD em

W . Para testar a exatidão

de algoritmos mistos, e também em estudos de transporte de elétrons de baixa energia

(

keV100 E < ), pode ser conveniente realizar simulações estritamente detalhadas. Com

este propósito, o PENELOPE permite ‘desligar’ a emissão de fótons de bremsstrahlung

suaves com energia menor do que 10 eV. Esta opção é ativada quando seleciona-se

W (caso em que o programa define eV10 W

) desconsiderando eventos de emissão

de bremsstrahlung suave e simulando eventos catastróficos (com

eV10 W > )

detalhadamente. A geração da deflexão angular em eventos artificiais é descontínua

quando a simulação de espalhamento elástico e inelástico se torna detalhada (isto é, quando

λλ

)(

e 0

W ).

O comprimento dos passos gerados pelo PENELOPE é sempre menor do que

max

s , limite

superior estabelecido pelo usuário. O código de simulação limita o tamanho do passo

colocando interações delta ao longo da trajetória da partícula. Estas são interações fictícias

que não alteram o estado da partícula, seu único efeito é interromper a seqüência de

operações da simulação, que requer a alteração dos valores das variáveis de controle

interno para permitir a retomada da simulação de uma maneira consistente. O uso de

passos de comprimento limitado é necessário devido à dependência energética das SCD’s

para interações suaves, porém esta não é a única razão para se limitar o tamanho do passo.

Uma vez que os valores das perdas de energia e deflexões nas dobradiças são obtidos a

partir de distribuições artificiais, o número de dobradiças no trajeto de uma partícula

primária deve ser estatisticamente suficiente, isto é, maior que ~10, para eliminar as

distorções introduzidas pela adoção de distribuições artificiais. Por essa razão, quando a

partícula atravessa um material estreito, é aconselhável o uso de um valor menor para

max

para garantir que o número de dobradiças que ocorrem dentro do meio seja suficiente. Para

garantir uma consistência interna,

max

s deve ser menor do que

)(

. Quando o valor

escolhido pelo usuário é maior, o código define

)(

max

= ; neste caso, aproximadamente

5% dos passos da amostra têm comprimento maior do que

max

s e são interrompidos por

uma interação delta. Isto torna a simulação um pouco lenta (~ 5%), mas garante que a

dependência energética de

)(h

seja considerada corretamente. No PENELOPE, o

parâmetro

max

s pode variar livremente durante o curso da simulação de uma mesma

trajetória, com base nisto, ao invés de usar o valor

max

s fornecido pelo usuário, o

PENELOPE usa um valor aleatório [obtido de uma distribuição triangular no intervalo

),0(

max

s ] cuja média é a metade do valor fornecido pelo usuário; isto é feito para eliminar

um artefato na distribuição de dose em profundidade de feixes paralelos de

elétrons/pósitrons próximo à superfície de entrada.

O estado de uma partícula, imediatamente após um evento, é definido pela sua posição

energia E e pelos cossenos diretores da sua direção de movimento

no referencial do

laboratório. Quando a energia da partícula torna-se menor do que o valor pré-definido

abs

E , assume-se que ela é localmente absorvida; considera-se que pósitrons aniquilam após

a absorção. A geração de trajetórias aleatórias para elétrons e pósitrons em meio material,

que pode ser constituído de várias regiões homogêneas de composições diferentes

separadas por superfícies (interfaces) bem definidas, é feita executando os seguintes

passos:

(i) Definir a posição inicial

r, energia cinética

e direção de movimento d

partícula primária.

(ii) Determinar a perda de energia ‘suave’ máxima

max

w ao longo de um passo e

definir o valor do inverso do livre caminho médio para eventos catastróficos.

Estes resultados dependem do valor adotado para

max

s , que pode variar ao longo

da trajetória.

(iii) Sortear a distância

s a ser percorrida até o próximo evento catastrófico (ou

interação delta) usando

max,

−

. ( 4.8 )

max

ss > , truncar o passo fazendo

max

(iv) Gerar o comprimento do passo

até a próxima dobradiça. Deixar a partícula

percorrer esta distância na direção

drrd

+← .

(v) Se a trajetória cruzar uma interface:

Parar a partícula no ponto de cruzamento (isto é, definir

r como igual à posição

deste ponto e fazer

igual à distância percorrida).

Ir para o passo (ii) para continuar a simulação no novo material, ou ir para (xi) se

o novo material for vácuo.

(vi) Simular a perda de energia e a deflexão na dobradiça. Este passo consiste de duas

ações:

a. Sortear o ângulo de deflexão polar

(

)

2cos1

−

da distribuição,

correspondente considerando a energia atual

. Sortear o ângulo de

espalhamento azimutal usando

πξ

. Realizar a rotação ),(

R do vetor

de acordo com as deflexões angulares polar e azimutal sorteadas (como

descrito na seção 3.2.2) para obter a nova direção

),(

φθ

R← .

b. Sortear a perda de energia

w devido às interações suaves ao longo do passo

s a partir da distribuição de probabilidade correspondente, e subtrair da

energia cinética

wEE

−

←

Estas duas ações são executadas em ordem aleatória para considerar a

dependência energética do ‘livre caminho médio suave’.

abs

EE < , ir para o passo (xi).

(vii) Deixar a partícula caminhar a distância

−

s na direção drrd

)(:

−+← .

(viii) Fazer igual ao passo (v).

(ix) Se no passo (iii) o tamanho do passo foi truncado, isto é,

max

ss = , simular uma

interação delta.

Ir para o passo (ii).

(x) Simular um evento catastrófico:

Sortear o tipo de interação de acordo com as probabilidades pontuais

,,,,

max,

)(

max,

)(

max,

)(

max,

)(

σσσσ

max,

)(

no caso de pósitrons. ( 4.9 )

Se o evento for uma interação delta, voltar para (ii).

Sortear o ângulo de espalhamento polar

e a perda de energia W da FDP

correspondente. Gerar o ângulo de espalhamento azimutal usando

πξ

Aplicar a rotação ),(

R no vetor d

para obter a nova direção: dd

),(

φθ

R← .

Reduzir a energia cinética da partícula:

WEE

−

←

Se, como resultado da interação, uma partícula secundária for emitida na direção

, com energia

abss

EE > , armazenar o estado inicial )

,,(

E dr .

Ir para (ii) se

abs

EE > .

(xi) Simular a trajetória dos elétrons secundários e fótons produzidos pela partícula

primária (ou por outras secundárias simuladas previamente) antes de iniciar uma

nova trajetória primária.

(SALVAT, 2001).

4.7 Simulação do chuveiro

Os chuveiros de elétrons e fótons são simulados chamando-se uma série de sub-rotinas que

serão brevemente descritas a seguir. A simulação é inicializada com a seguinte sub-rotina:

Sub-rotina PEINIT (EPMAX, NMAT, IRD, IWR, INFO). A subrotina PEINIT lê os

arquivos de dados dos diferentes materiais, avalia as propriedades de espalhamento

relevantes e prepara tabelas para consultas relativas às grandezas dependentes da energia

que são usadas durante a simulação. Seus argumentos de entrada são:

EPMAX Energia máxima das partículas simuladas. Note que se as partículas

primárias forem pósitrons com energia cinética inicial EP, a energia

máxima dos fótons produzidos por aniquilação pode ser próxima de (mas

menor que)

(

)

EP21.1EPMAX cm

+= ; neste caso em especial, a energia

máxima é maior do que a energia cinética inicial.

NMAT Número de materiais envolvidos (menor ou igual à MAXMAT).

IRD Unidade de entrada.

IWR Unidade de saída.

INFO Determina a quantidade de informação que é escrita no arquivo de saída.

É mínima para

0INFO

, aumentando progressivamente o nível de

detalhamento para ,2 ,1INFO

etc.

Sub-rotina CLEANS. Origina a pilha secundária.

Sub-rotina START. Para elétrons e pósitrons, esta sub-rotina força o evento de interação

seguinte a ser um evento artificial suave. Ela deve ser chamada antes de começar uma nova

trajetória, primária ou secundária, e também quando uma trajetória cruza uma interface.

Chamar START é estritamente necessário para elétrons e pósitrons; para fótons esta sub-

rotina não tem efeito físico. Porém, é aconselhável chamar START para qualquer tipo de

partícula uma vez que nela é checado se a energia da partícula está dentro do limite

esperado, ajudando assim a detectar erros no programa.

Sub-rotina JUMP (DSMAX, DS). Determina o comprimento DS do segmento até o

próximo evento de interação.

O parâmetro de entrada DSMAX define o comprimento máximo do passo para

elétrons/pósitrons; para fótons ele não tem efeito. Como mencionado anteriormente, para

limitar o tamanho do passo, o PENELOPE insere interações delta ao longo da trajetória da

partícula.

Sub-rotina KNOCK (DE, ICOL). Simula um evento de interação, determina a nova energia

e direção de movimento, e armazena o estado inicial das partículas secundárias, caso sejam

geradas. Os argumentos de saída são:

DE energia depositada no curso do evento.

ICOL tipo de evento que foi simulado, de acordo com a seguinte convenção,

• Elétrons ( 1KPA

= )

ICOL = 1, evento artificial suave (dobradiça).

= 2, colisão elástica catastrófica.

= 3, colisão inelástica catastrófica.

= 4, emissão de bremsstrahlung catastrófica.

• Fótons ( 2KPA

= )

ICOL = 1, espalhamento coerente (Rayleigh).

= 2, espalhamento incoerente (Compton).

= 3, absorção fotoelétrica.

= 4, produção de pares.

• Pósitrons (

3KPAR = )

ICOL = 1, evento artificial suave (dobradiça).

= 2, colisão elástica catastrófica.

= 3, colisão inelástica catastrófica.

= 4, emissão de bremsstrahlung catastrófica.

= 5, aniquilação.

onde o parâmetro KPAR indica o tipo de partícula. Para elétrons e pósitrons

7ICOL

corresponde a interações delta. O valor

6ICOL

é usado para interações auxiliares (um

mecanismo de interação adicional que pode ser definido pelo usuário para simular, por

exemplo, interações fotonucleares).

Sub-rotina SECPAR (LEFT). Fornece o estado inicial de uma partícula secundária e a

remove da pilha secundária. O valor de saída LEFT é o número de partículas secundárias

que permanecem na pilha.

Sub-rotina STORES (E, X, Y, Z, U, V, W, WGHT, KPAR, ILB). Armazena uma partícula

na pilha secundária. Os argumentos são as variáveis de estado da partícula:

E energia da partícula (energia cinética para elétrons e pósitrons).

X, Y, Z coordenadas da posição.

U, V, W cossenos diretores da direção de movimento.

WGHT em simulações análogas, esta variável não tem efeito. Ao usar métodos de

redução de variância, pode ser usada para armazenar o peso da partícula.

KPAR tipo de partícula (1: elétron, 2: fóton, 3: pósitron).

ILB(5) arranjo auxiliar de 5 rótulos que descrevem a origem das partículas.

A seqüência de sub-rotinas chamadas para gerar uma trajetória aleatória independe do tipo

de partícula que está sendo simulada. A geração dos chuveiros aleatórios procede como a

seguir (FIGURA 4.2):

(i) Determinar o estado inicial da partícula primária, isto é, define-se os valores das

variáveis de estado KPAR, E, Z)Y,X,(

r e W)V,U,(

=d . Especificar o ‘corpo’

e o material em que a partícula move, definindo-se os valores de IBODY e MAT,

respectivamente. Opcionalmente, define-se os valores de WGHT e 5):ILB(1 .

(ii) Chamar a sub-rotina CLEANS para inicializar a pilha secundária.

(iii) Chamar a sub-rotina START para iniciar a simulação da trajetória.

(iv) Chamar a sub-rotina JUMP (DSMAX, DS) para determinar o comprimento DS do

próximo segmento de trajetória (para elétrons e pósitrons, DS não pode exceder o

valor de entrada DSMAX).

(v) Calcular a posição do evento seguinte.

• Se a trajetória cruza uma interface, ela é interrompida na posição em que a

interseção ocorre, reduzindo-se DS.

Mudar para o novo material redefinindo as variáveis IBODY e MAT.

Quando a partícula escapa do sistema, a simulação é finalizada; atualizar os

contadores e ir para o passo (vii).

Ir para o passo (iii).

(vi) Chamar a sub-rotina KNOCK (DE, ICOL) para simular o próximo evento.

• Se a energia tornar-se menor do que EABS (KPAR, MAT), finalizar a

trajetória, atualizar os contadores e ir para o passo (vii).

• Ir para o passo (iv).

(vii) Chamar a sub-rotina SECPAR (LEFT) para iniciar a trajetória de uma partícula da

pilha secundária (esta partícula é então automaticamente removida da pilha).

• Se LEFT > 0, ir para o passo (iii). Já tendo sido determinado o estado inicial

da partícula secundária.

• Se LEFT = 0, a simulação do chuveiro produzido pela partícula primária terá

terminado. Ir para o passo (i) para gerar uma nova partícula primária (ou

terminar a simulação se já tiver sido simulado um número suficiente de

chuveiros).

Note que as sub-rotinas JUMP e KNOCK mantêm as coordenadas de posição inalteradas;

as posições dos sucessivos eventos são deduzidas a partir do deslocamento de

comprimento DS ao longo da direção de movimento após cada chamada da sub-rotina

JUMP. A energia da partícula é automaticamente reduzida pela sub-rotina KNOCK, que

também modifica a direção de movimento de acordo com os ângulos de espalhamento do

evento simulado. Assim, na saída de KNOCK, os valores da energia

, posição

Z)Y,X,(=

r e direção de movimento W)V,U,(

=d definem o estado da partícula

imediatamente após o evento de interação (SALVAT, 2001).

4.8 Estabilidade do Algoritmo de Simulação

No caso da simulação do transporte de elétrons/pósitrons de baixa energia (da ordem de

500 keV ou menos), os parâmetros relevantes são

abs

E ,

C ,

W e

max

s , uma vez que

não é efetivo e a emissão radiativa é desprezível (eventos de bremsstrahlung catastróficos

ocorrem muito raramente e, então,

W não tem nenhuma influência). O valor do

parâmetro

max

s é importante para garantir a confiabilidade dos resultados; uma maneira

segura de se definir este parâmetro é fazê-lo igual a um décimo do valor “esperado” para o

comprimento total da trajetória ou menos. Uma vez que os valores de

abs

E e

W são

ditados pelas características do experimento considerado, o único parâmetro crítico, que

influencia diretamente a velocidade da simulação, é

C . Como mencionado anteriormente,

o PENELOPE aceita valores de

C entre 0 (simulação detalhada de espalhamentos

elásticos) e 0,2. Na prática, o valor de

C não influencia a precisão dos resultados quando

aos outros parâmetros são dados valores “seguros” (SALVAT, 2001). Neste trabalho foram

utilizados os valores indicados como seguros no manual do PENELOPE (SALVAT, 1996

citado por STEWART, 2002): C

= 0,001 e C

= 0,01. Para as energias de corte foi

escolhido o menor valor com que o PENELOPE trabalha: W

= W

= 100 eV.

Inicializa o PENELOPE

Inicia um novo chuveiro

Estado inicial

sim

A trajetória cruza uma interface?

Modificar DS para terminar

rimeira interface

Mudar para o novo corpo

Saiu?

Figura 4.2: Fluxograma do programa para simulação de chuveiros de elétrons e

fótons com o PENELOPE.

5 O cálculo dos alcances

Para o cálculo dos alcances foram desenvolvidos algoritmos na linguagem de programação

FORTRAN que acessam os arquivos de dados em que foram gravadas as informações a

respeito de toda a história de cada partícula simulada para obter as informações necessárias

para estes cálculos. A título de ilustração, encontra-se abaixo um dos algoritmos

desenvolvidos para este fim.

C*******************************************************************

C SUB-ROTINA SOMA

C*******************************************************************

C Calcula a média do comprimento das trajetórias das NTOT partículas

C simuladas com uma mesma energia inicial E0 (partículas primárias).

SUBROUTINE SOMA(R)

IMPLICIT DOUBLE PRECISION (A-H,O-Z), INTEGER*4 (I-N)

DIMENSION ILB(5)

OPEN(60,FILE='temp.dat',FORM='unformatted',STATUS='scratch')

READ(65) E0,NTOT

10 CONTINUE

READ(65,END=30) E,DE,X,Y,Z,U,V,W,DS,KPAR,ICOL,ILB(1),ILB(2),

1 ILB(3),ILB(4),ILB(5)

Ri=0

IF (E.EQ.E0) THEN

Ri=Ri+DS

ELSE

GO TO 10

END IF

20 CONTINUE

READ(65,END=30) E,DE,X,Y,Z,U,V,W,DS,KPAR,ICOL,ILB(1),ILB(2),

1 ILB(3),ILB(4),ILB(5)

IF (E.NE.0) THEN

Ri=Ri+DS

GO TO 20

END IF

Ri=Ri+DS

WRITE(60) Ri

GO TO 10

30 CONTINUE

REWIND(60)

SOM=0

40 CONTINUE

READ(60,END=50) Ri

SOM=SOM+Ri

GO TO 40

50 CONTINUE

R=SOM/NTOT

CLOSE(60)

RETURN

END

A sub-rotina SOMA calcula o valor médio, R, dos comprimentos das trajetórias de NTOT

elétrons simulados com uma mesma energia inicial E0. Isto é feito somando, para cada

elétron primário, todos os deslocamentos DS entre duas interações consecutivas desde

E = E0 até E = 0, ou seja, do ponto inicial até o ponto final da trajetória deste elétron,

gravando o valor desta soma (que é igual ao comprimento da trajetória) em um segundo

arquivo para cada um dos NTOT elétrons primários da simulação. Este arquivo é então lido

e o comprimento médio das trajetórias destes NTOT elétrons calculado. Esta média por sua

vez é gravada em um terceiro arquivo onde são registrados os valores de R das dez

simulações realizadas. Daí é calculado o alcance R propriamente dito, com o desvio padrão

(incerteza) correspondente: a partir da média destes dez valores obtidos em dez simulações

com NTOT chuveiros. Todos os demais alcances foram calculados desta mesma maneira,

fazendo-se a média de dez valores obtidos em simulações com um mesmo número NTOT

de chuveiros.

Com as constantes interações que sofrem ao longo de sua trajetória, os elétrons penetram

no meio material até que sua energia cinética entre em equilíbrio térmico com as partículas

do meio, estabelecendo um “alcance” no meio absorvedor. Devido às dificuldades

experimentais na determinação do “alcance” de elétrons de baixa energia, pode-se definir

diferentes alcances utilizando variações da definição. Uma partícula carregada, como o

elétron, percorre um caminho aproximadamente bem definido no material antes de dissipar

toda sua energia. Devido à natureza estatística do processo de perda de energia, partículas

com uma mesma energia inicial terão trajetórias de comprimento e forma diferentes. Por

este motivo, os alcances calculados aqui são todos uma média dos valores determinados

para um grande número de elétrons primários com uma mesma energia inicial. Estes

alcances médios, incluindo o alcance R, foram calculados para cada energia considerada –

dentro do intervalo de 110 eV e 50 keV – baseados nas seguintes definições:

Alcance R:

É a soma de todos os deslocamentos do elétron primário entre duas interações consecutivas

desde o ponto inicial até o ponto final da trajetória deste, ou seja, o alcance R é igual ao

comprimento da trajetória do elétron (FIGURA 4.3):

FIGURA 4.3: Representação do alcance R.

∑

−=

R rr

Este alcance foi determinado como a seguir:

(i) Calculou-se o tamanho da trajetória dos NTOT elétrons primários de uma

simulação, somando todos os deslocamentos DS do elétron desde E = E0 até E = 0;

(ii) Calculou-se então o comprimento médio das trajetórias destes NTOT elétrons;

(iii) Repetiu-se os passos (i) e (ii) dez vezes e calculou-se a média e o desvio padrão

destes dez valores.

Alcance R

Levando em conta que a fonte de elétrons considerada neste trabalho é isotrópica, o

alcance R

foi definido como sendo igual ao raio do volume dentro do qual 90% da

energia inicial do elétron primário foi depositada. Foi determinado como a seguir:

(i) Em toda simulação de NTOT chuveiros, o código de Monte Carlo PENELOPE

calcula uma distribuição de dose em profundidade média por elétron. A partir desta

distribuição de dose em profundidade determinou-se o raio em que 90% da energia

do elétron primário foi depositada;

(ii) Repetiu-se o passo (i) dez vezes e calculou-se a média e o desvio padrão destes dez

valores.

Alcance R

Analogamente ao alcance R

, o alcance R

indica o raio do volume dentro do qual 95%

da energia inicial do elétron primário foi depositada. Seu cálculo foi similar ao do

alcance R

Alcance projetado, R

proj

É a penetração máxima do elétron primário em uma determinada direção, em geral a

direção de incidência (FIGURA 4.4). No caso deste trabalho, em que a fonte de elétrons é

isotrópica, R

proj

não depende da direção considerada. Este alcance foi determinado como a

seguir:

FIGURA 4.4: Representação do alcance projetado, R

proj

(i) Procurou-se o máximo, no arquivo de dados, de uma das coordenadas (x, y ou z)

para cada um dos NTOT elétrons primários de uma simulação, ou seja, determinou-

se a penetração máxima dos elétrons primários em uma determinada direção;

(ii) Calculou-se o valor médio destes máximos;

(iii) Repetiu-se os passos (i) e (ii) dez vezes e calculou-se a média e o desvio padrão

destes dez valores.

Deslocamento total, DT:

É a distância entre os pontos inicial (r

) e final (r

) da trajetória da partícula (FIGURA 4.5),

isto é,

DT −=

O deslocamento total do elétron foi calculado como a seguir:

(i) Determinou-se as coordenadas dos pontos em que E = E0 e E = 0, ou seja, dos

pontos inicial e final da trajetória dos NTOT elétrons primários de uma simulação;

(ii) Calculou-se a distância entre estes pontos para os NTOT elétrons;

(iii) Repetiu-se os passos (i) e (ii) dez vezes e calculou-se a média e o desvio padrão

destes dez valores.

Alcance ponderado, R

pond

É a soma dos deslocamentos DS, desde o ponto inicial (E = E0) até o ponto final (E = 0) da

trajetória do elétron primário, ponderada pela energia depositada em cada um destes

deslocamentos, isto é,

∑

∆

∆−

iii

pond

onde

r é a posição e

E∆ a energia depositada no i-ésimo evento. Este alcance foi

determinado como a seguir:

FIGURA 4.5: Representação do deslocamento total, DT.

(i) Fez-se a soma dos deslocamentos DS, do início ao fim da trajetória do elétron,

multiplicados pela energia depositada ao longo de cada um destes passos. Ao

mesmo tempo acumulou-se em um contador estes valores de energia depositada;

(ii) Calculou-se então para cada um dos NTOT elétrons a razão entre estas somas;

(iii) Finalmente, calculou-se a média dos alcances obtidos no item anterior;

(iv) Repetiu-se os passos (i), (ii) e (iii) dez vezes e calculou-se a média e o desvio

padrão destes dez valores.

Os valores obtidos neste trabalho foram comparados com valores experimentais e teóricos

determinados em tecido ou em água. A comparação com valores de alcance em água é

justificável porque a densidade do tecido utilizado nas simulações (tecido mole ICRP) é

igual a um.

Conhecer o alcance de elétrons de baixa energia é relevante, dentre outros motivos, para se

determinar a extensão sobre a qual, elétrons colocados em movimento no interior de um

meio biológico, podem produzir danos relevantes às entidades biológicas irradiadas. Como

ocorre, por exemplo, no uso de anticorpos monoclonais para o tratamento dos mais

variados tipos de neoplasias. Além disto, dispondo de uma variedade de definições para o

alcance de um elétron pode-se determinar qual a definição mais adequada a uma

determinada aplicação. Por exemplo, o alcance R95 poderia ser utilizado quando o

interesse está na região de deposição da energia inicial do elétron, uma vez que este

alcance é uma medida do raio dentro do qual 95% da energia inicial do elétron é

depositada.

Na próxima seção serão apresentados os resultados obtidos neste trabalho para o cálculo de

alcances de elétrons com energias iniciais de 110, 200, 300, 400, 500, 600, 700, 800 e

900 eV, e 1, 5, 10, 25 e 50 keV provenientes de uma fonte monoenergética isotrópica

localizada no interior de um meio tecido equivalente homogêneo semi-infinito.

100

RESULTADOS

As FIGURAS 5.1 e 5.2 apresentam as curvas dos alcances em tecido biológico calculados

neste trabalho versus a energia inicial do elétron. Para uma melhor visualização, a região

de mais baixa energia (de 110 eV a 1 keV) foi plotada separadamente (FIGURA 5.2).

Todos valores calculados têm incerteza máxima de 1% e estão compreendidos no intervalo

em torno de aproximadamente 5 Å e 40 µm.

Como era esperado, o alcance R, dado pelo comprimento total da trajetória do elétron, é

maior do que todos os outros alcances. Já o alcance projetado, R

proj

, dado pela maior

profundidade de penetração do elétron primário em uma determinada direção, é o menor

dentre eles. Curiosamente, os alcances DT e R

apresentam valores comparáveis,

principalmente para energias entre 110 eV e 1 keV. O alcance ponderado em energia,

pond

R , foi o único que apresentou comportamento peculiar em função da energia inicial do

0 1020304050

0,0

1,0x10

-3

2,0x10

-3

3,0x10

-3

4,0x10

-3

5,0x10

-3

Alcance (cm)

Energia (keV)

Projetado

R90

R95

FIGURA 5.1: Alcances de elétrons com energias entre 110 eV e

50 keV em tecido biológico.

proj

Energia inicial do elétron (keV)

101

elétron, por este motivo ele não foi incluído nas FIGURAS 5.1 e 5.2. Como pode ser

observado na FIGURA 5.3, diferentemente dos outros alcances, o alcance ponderado

aumenta linearmente com o aumento da energia inicial do elétron e, além disto, seus

valores são muito pequenos se comparados com os demais alcances – sendo o maior valor,

para energia inicial de 50 keV, inferior a 0,1 µ.

0 200 400 600 800 1000

0,0

5,0x10

-7

1,0x10

-6

1,5x10

-6

2,0x10

-6

2,5x10

-6

3,0x10

-6

3,5x10

-6

4,0x10

-6

4,5x10

-6

5,0x10

-6

5,5x10

-6

Alcance (cm)

Energia (eV)

Projetado

R90

R95

FIGURA 5.2: Alcances de elétrons com energias entre 110 eV e 1 keV em tecido biológico.

Energia inicial do elétron (eV)

0 1020304050

0,0

1,0x10

-6

2,0x10

-6

3,0x10

-6

4,0x10

-6

5,0x10

-6

6,0x10

-6

Alcance (cm)

Energia (keV)

Alcance Ponderado

FIGURA 5.3: Alcance ponderado versus energia inicial do

elétron para elétrons com energias entre 110 eV e 50 keV.

102

A regressão linear da curva apresentada na FIGURA 5.3 fornece a seguinte expressão para

o alcance ponderado em centímetros:

pond

108

1007607,11002088,9

−−

×+×=

onde E é a energia inicial do elétron em elétron-volts. Dividindo o alcance ponderado

∑∑

∆∆−=

+ iiipond

EER rr

pelo alcance

∑

−=

+ ii

R rr

, temos:

total

iii

pond

∆

−∆−

∑

++ i1i1

rr rr

Isto é, a razão RR

pond

/ é igual ao valor médio da fração da energia do elétron depositada

por passo. O gráfico de RR

pond

/ versus a energia inicial do elétron (FIGURA 5.4), mostra

que esta razão decresce com o aumento da energia, isto é, quanto maior for a energia do

elétron, menor será a fração desta energia depositada em cada passo. Isto indica que é no

final da trajetória do elétron, onde sua energia e seu livre caminho médio são menores, que

a maior parte da energia inicial do elétron é depositada.

100 1000 10000

-5

pond

/R (%)

log(E) (keV)

(Alcance Ponderado / Alcance R)

FIGURA 5.4: Razão entre os alcances ponderado, R

pond

, e

R em função da energia inicial do elétron.

Energia inicial do elétron (eV)

103

As FIGURAS 5.5 e 5.6 comparam o alcance R com alguns alcances encontrados na

literatura determinados teórica e/ou experimentalmente: o alcance determinado por Cole,

Cole

, por Lea, R

Lea

, e os alcances CSDA1

e CSDA2. O alcance R

Cole

é um alcance semi-

empírico, já os demais alcances (R

Lea

, CSDA1 e CSDA2) são todos alcances teóricos

calculados utilizando a aproximação de frenamento contínuo. Como pode ser observado

nestas figuras, o alcance R apresenta valores parecidos com os destes alcances, exceto para

energias menores do que 500 eV em que esta diferença aumenta consideravelmente,

principalmente em relação ao alcance R

Cole

Cole (1969) determinou alcances para elétrons com energias entre 20 eV e 50 keV

medindo a transmissão de partículas através de lâminas de colódio. Este alcance foi

definido como sendo igual à espessura em que 5% das partículas são transmitidas através

do material absorvedor (isto é, 95% das partículas incidentes são absorvidas em um

material de espessura igual a este alcance). A partir dos alcances medidos, Cole

determinou a seguinte relação empírica:

007.0)367.0(0431.0

77.1

−+= ER

Cole

onde E é dado em keV e o alcance em

g/cm 100

, que corresponde a 1

para

densidade unitária. A estes valores está associada uma incerteza de aproximadamente 10%.

CSDA = Continuous-Slowing-Down-Approximation.

0 1020304050

0,0

1,0x10

-3

2,0x10

-3

3,0x10

-3

4,0x10

-3

5,0x10

-3

Alcance (cm)

Energia (keV)

Cole

CSDA1

CSDA2

FIGURA 5.5: Com

ara

ão do alcance

com os alcances

Cole

CSDA1 e CSDA2.

Cole

Energia inicial do elétron (keV)

104

Energia

(eV)

R (cm)

Cole

(cm)

CSDA1

(cm)

CSDA2

(cm)

Lea (cm)

Diferença

% de

Cole

Diferença

% de

CSDA1

Diferença

% de

CSDA2

Diferença

% de Lea

110 1,64E-07 4,20E-07 5,00E-07 -- 4,55E-07 61 67 -- 64

200 4,01E-07 8,79E-07 8,00E-07 -- 7,81E-07 54 50 -- 49

300 7,61E-07 1,40E-06 1,15E-06 -- 1,19E-06 46 34 -- 36

400 1,19E-06 2,00E-06 1,55E-06 -- 1,66E-06 40 23 -- 28

500 1,70E-06 2,65E-06 2,00E-06 -- 1,98E-06 36 15 -- 14

600 2,27E-06 3,36E-06 2,51E-06 -- 2,55E-06 32 9 -- 11

700 2,91E-06 4,13E-06 3,06E-06 -- 3,10E-06 30 5 -- 6

800 3,61E-06 4,96E-06 3,66E-06 -- 3,80E-06 27 1 -- 5

900 4,37E-06 5,85E-06 4,30E-06 -- 4,52E-06 25 2 -- 3

1000 5,19E-06 6,80E-06 5,00E-06 -- 5,30E-06 24 4 -- 2

5000 7,71E-05 8,37E-05 8,00E-05 -- 7,99E-05 8 4 -- 4

10000 2,57E-04 2,70E-04 2,87E-04 2,54E-04 2,91E-04 5 10 1,4 12

25000 1,29E-03 1,32E-03 1,23E-03 1,28E-03 1,71E-03 2 5 0,8 25

50000 4,39E-03 4,44E-03 4,25E-03 4,36E-03 6,72E-03 1 3 0,7 35

TABELA 5.1: Alcances R, R

Cole

, R

Lea

, CSDA1 e CSDA2 e diferença percentual de R e

relação a estes alcances. Os valores de R

Lea

relativos a energias maiores que 1 keV fora

obtidos por extrapolação.

Cole

Lea

(cm)

FIGURA 5.6: Destaque da região de baixa energia (110 eV a 1 keV). Comparação do

alcance R com os alcances R

Cole

, CSDA1 e R

Lea

0 200 400 600 800 1000

0,0

1,0x10

-6

2,0x10

-6

3,0x10

-6

4,0x10

-6

5,0x10

-6

6,0x10

-6

7,0x10

-6

Alcance (cm)

Energia (eV)

Cole

CSDA1

Lea

Energia inicial do elétron (eV)

Cole

105

Os alcances CSDA1 foram calculados a partir de fórmulas analíticas para o poder de

frenamento total de elétrons em água (TURNER, 1995), com extrapolação para baixas

energias. Já os alcances CSDA2 foram retirados da base de dados ESTAR

(BERGER, 2000) que calcula tabelas com valores de poder de frenamento e alcance para

elétrons de acordo com os métodos

descritos no ICRU 37 (ICRU, 1984). Neste método, o

poder de frenamento de colisão é calculado a partir da teoria de Bethe (1930, 1932) (citado

por BERGER, 2000), com correção para o efeito da densidade de acordo com a calculada

por Sternheimer (1952, 1982) (citado por BERGER, 2000), e o poder de frenamento

radiativo, a partir da combinação das seções de choque de bremsstrahlung descritas por

Seltzer e Berger (1985) (citado por BERGER, 2000) com fórmulas analíticas e resultados

numéricos de Pratt et al (1977) (citado por BERGER, 2000). As incertezas associadas a

estes valores, que foram determinados para tecido biológico, estão em torno de 3%. Como

os alcances CSDA1 e CSDA2 foram determinados em água e tecido respectivamente, eles

apresentam uma pequena diferença que varia entre aproximadamente 3 e 12%.

Para baixas energias, além dos resultados de Cole, R

Cole

, e da extrapolação CSDA1, estão

representados no gráfico (FIGURA 5.6) os alcances calculados por Lea (1952), R

Lea

, a

partir da teoria de Bethe para perdas de energia, assumindo a validade da extrapolação da

expressão até baixas energias. Estes cálculos, que foram feitos para água, estão corrigidos

para absorção em tecido. Vale salientar que, enquanto a aproximação de frenamento

contínuo (CSDA) é fisicamente apropriada para partículas carregadas pesadas, ela não é

sempre realista para elétrons, que podem perder grandes frações de sua energia em uma

única colisão (TURNER, 1995). Além disto, os valores dos alcances calculados usando a

equação de Bethe abaixo da região de aplicabilidade tornam-se extremamente dependentes

da seleção do valor do potencial de excitação/ionização médio, I (COLE, 1969). A difícil

definição do potencial de excitação/ionização se deve à não existência de correções para as

camadas atômicas que são necessárias quando a velocidade do elétron incidente é pequena

se comparada às velocidades orbitais dos elétrons atômicos. Em função desta limitação,

estes cálculos que utilizam a aproximação de frenamento contínuo são, usualmente,

restritos a energias maiores que 10 keV (BERGER, 2000).

Estes métodos foram desenvolvidos por um comitê formado pela International Commission on Radiation

Units and Measurements. Eram membros deste comitê H. H. Anderson, M. J. Berger, H. Bichsel,

J. A. Dennis, M. Inokuti, D. Powers, S. M. Seltzer, D. Thwaites, J. E. Turner e D. E. Watt.

106

Na aproximação de frenamento contínuo, assume-se que a partícula perde energia

continuamente ao longo de sua trajetória independentemente das deflexões sofridas no

percurso (isto é, independente da trajetória ser retilínea ou não) a uma taxa que depende da

energia. Sendo assim, os alcances calculados utilizando a aproximação de frenamento

contínuo (CSDA1, CSDA2 e R

Lea

) dão uma medida aproximada do comprimento da

trajetória da partícula em função da energia inicial desta. Como esperado, os valores

obtidos para o alcance R é bastante próximo aos valores dos alcances CSDA1, CSDA2 e

Lea

e à medida que a energia diminui, a discrepância entre os valores aumenta,

possivelmente em conseqüência da restrição do uso da aproximação de frenamento

contínuo no cálculo de alcances de elétrons com energia inicial menor que 10 keV. Na

região de mais altas energias (de 1 a 50 keV), o alcance R difere de 1 a 8% do R

Cole

, de 3 a

10% do alcance CSDA1 e de 0,7 a 1,4% do CSDA2 (TABELA 5.1). Em baixas energias

(de 110 eV a 1 keV) as diferenças aumentam consideravelmente à medida que a energia

diminui, como pode ser visto na TABELA 5.1, chegando a 61% em relação ao R

Cole

, 67%

em relação ao CSDA1 e a 64% em relação ao R

Lea

. Nota-se que, além de serem grandes as

diferenças do alcance R em relação aos demais em energias mais baixas, estas diferenças

são da mesma ordem (~ 60%) em 110 e 200 eV, mesmo sendo o alcance R

Cole

determinado

experimentalmente e não teoricamente como os alcances CSDA1 e R

Lea

. Isto pode ter

origem no fato de que, embora o PENELOPE possa transportar elétrons até energias

comparáveis a 100 eV, as seções de choque e a física de maneira geral, para energias

abaixo de 500 eV têm incertezas maiores do que para altas energias, pois nesta região não

são considerados os refinamentos nas seções de choque diferenciais. O alcance R é

portanto uma boa estimativa para o tamanho da trajetória, em tecido, para elétrons com

energia inicial maior do que 500 eV.

A FIGURA 5.7 apresenta os resultados dos alcances R

, R

, DT e R

proj

comparados com

os dados experimentais

determinados por Davis (1954), R

Davis

. Davis mediu os alcances de

elétrons com energias entre 600 eV e 2 keV utilizando uma enzima biologicamente ativa

como “detetor”: quando uma camada uniforme desta enzima é irradiada, as moléculas

desta enzima situadas dentro do alcance dos elétrons são inativadas; não sendo afetadas

pela radiação as enzimas localizadas além deste alcance. Após a irradiação de um volume

conhecido desta enzima com elétrons monoenergéticos, era determinada a massa

Na realidade, os alcances estabelecidos por Davis foram obtidos ajustando-se uma curva aos dados das

medidas experimentais. Portanto, estes resultados são semi-empíricos.

107

correspondente à quantidade de enzimas inativadas. Uma vez conhecidas a área irradiada, a

densidade e a massa da fração de enzimas inativas, determina-se então o alcance. Os

valores utilizados aqui foram corrigidos para densidade unitária simplesmente

multiplicando os resultados de Davis pela densidade de 1,3 g/cm

da enzima utilizada

. A

incerteza destes dados é de aproximadamente 26%.

Comparado com os valores de Davis na faixa de 200 eV a 1 keV (TABELA 4.2), o alcance

DT apresenta diferenças no intervalo de aproximadamente 0,3 a 10,5%, o R

de 1,4 a

22,1% e o R

de 0,1 a 24,5%. Para a energia de 110 eV, estas diferenças são maiores para

os alcances DT e R

, sendo iguais a 48 e 38% respectivamente. A diferença percentual

média destes alcances em relação a R

Davis

é de 10,1% para o alcance DT, de 14,9% para o

e de 12,5% para o R

, portanto, o alcance DT é o que mais se aproxima do alcance

experimental de Davis.

Ao realizar uma correção desta maneira, considera-se hipoteticamente que o comportamento dos eventos de

espalhamento é linear em relação à densidade do meio. Cumpre observar que esta prática, apesar de não ser a

mais adequada, é frequentemente utilizada neste tipo de comparação.

0 200 400 600 800 1000

0,0

5,0x10

-7

1,0x10

-6

1,5x10

-6

2,0x10

-6

2,5x10

-6

3,0x10

-6

3,5x10

-6

4,0x10

-6

Alcance (cm)

Energia (eV)

Projetado

R90

R95

Davis

proj

Davis

TABELA 4.2: Comparação dos alcances DT, R

proj

, R

e R

com o alcance de Davis, R

Davis

, com as respectivas diferenças percentuais

em relação a este alcance.

Energia

(eV)

DT (cm)

Projetado

(cm)

R90 (cm) R95 (cm) Davis (cm)

Diferença %

de DT

Diferença %

de R90

Diferença %

de R95

110 8,94E-08 4,68E-08 1,05E-07 1,33E-07 1,72E-07 48,0 38,9 23,0

200 2,48E-07 1,35E-07 2,50E-07 3,22E-07 2,72E-07 9,0 8,3 18,3

300 4,94E-07 2,74E-07 4,48E-07 5,65E-07 4,54E-07 8,8 1,4 24,5

400 7,85E-07 4,40E-07 6,74E-07 8,66E-07 7,10E-07 10,5 5,1 22,0

500 1,12E-06 6,33E-07 9,71E-07 1,18E-06 1,04E-06 7,6 6,7 13,6

600 1,50E-06 8,50E-07 1,30E-06 1,57E-06 1,45E-06 3,6 10,1 8,6

700 1,92E-06 1,09E-06 1,59E-06 1,92E-06 1,93E-06 0,3 17,4 0,1

800 2,38E-06 1,35E-06 1,99E-06 2,39E-06 2,46E-06 3,1 18,9 2,7

900 2,88E-06 1,64E-06 2,39E-06 2,87E-06 3,00E-06 4,1 20,4 4,6

1000 3,42E-06 1,95E-06 2,84E-06 3,37E-06 3,64E-06 6,0 22,1 7,3

5000 5,16E-05 2,93E-05 4,26E-05 5,06E-05 -- -- -- --

10000 1,74E-04 9,86E-05 1,42E-04 1,69E-04 -- -- -- --

25000 8,87E-04 4,99E-04 7,19E-04 8,37E-04 -- -- -- --

50000 3,03E-03 1,70E-03 2,44E-03 2,84E-03 -- -- -- --

proj

Davis

109

Observando o gráfico que mostra a razão entre os alcances DT e R em função da energia

inicial do elétron (em todo espectro considerado, de 110 eV a 50 keV) (FIGURA 4.8),

percebe-se que esta razão aumenta rapidamente de 110 eV até 1 keV. A partir de 1 keV, a

razão DT / R cresce mais lentamente e parece tender a se estabilizar quando o valor desta

se aproxima de 69%. Este comportamento indica que à medida que a energia aumenta, a

trajetória do elétron tende a ser menos tortuosa, fazendo com que diminua a diferença entre

a distância DT, entre pontos inicial e final da trajetória, e o tamano desta, R. No caso da

razão entre os alcances projetado, R

proj

, e R (FIGURA 4.9), nota-se o mesmo

comportamento: um crescimento rápido da razão entre 110 eV e 1 keV e um crescimento

mais lento a partir daí que parece tender a se estabilizar ao se aproximar de 39%. A

diferença entre os alcances R

proj

e R é aproximadamente constante e fica em torno de 63%,

exceto em 110 eV que é de ~ 72% (TABELA 5.2), o que concorda razoavelmente com a

estimativa de 60% de Williams (WILLIAMS, 1930 citado por DAVIS, 1954). Interessante

observar que, assim como o alcance DT, o alcance projetado tende a manter uma diferença

constante em relação ao alcance R para energias maiores que 1 keV e, neste caso, o

alcance projetado é praticamente igual à metade do tamanho da trajetória do elétron, R.

FIGURA 5.8: Razão entre os alcances DT e R em função da energia

inicial do elétron.

0 1020304050

DT/ R (%)

Energia (keV)

DT/ R

0,1 1 10

DT/ R (%)

log (E) (keV)

Energia inicial do elétron (keV)

110

Energia (eV) Projetado (cm)

R (cm) Projetado/R Diferença %

110 4,68E-08

1,64E-07 0,28506 72

200 1,35E-07

4,01E-07 0,33593 66

300 2,74E-07

7,61E-07 0,36072 64

400 4,40E-07

1,19E-06 0,36887 63

500 6,33E-07

1,70E-06 0,37231 63

600 8,50E-07

2,27E-06 0,3737 63

700 1,09E-06

2,91E-06 0,37435 63

800 1,35E-06

3,61E-06 0,37499 63

900 1,64E-06

4,37E-06 0,3753 63

1000 1,95E-06

5,19E-06 0,37539 63

5000 2,93E-05

7,71E-05 0,38039 62

10000 9,86E-05

2,57E-04 0,38339 62

25000 4,99E-04

1,29E-03 0,38523 62

50000 1,70E-03

4,39E-03 0,38776 61

TABELA 5.2: Comparação entre os alcances projetado R

proj

e R.

proj

(cm) R

proj

/ R

0 1020304050

proj

/ R (%)

Energia (keV)

Alcance projetado / Alcance R

0,1 1 10

proj

/ R (%)

log (E) (keV)

FIGURA 5.9: Razão entre os alcances projetado R

proj

e R em

função da energia inicial do elétron.

Energia inicial do elétron (keV)

111

CONCLUSÕES

O alcance R é uma boa estimativa para o tamanho da trajetória, em tecido, de elétrons com

energia inicial maior do que 500 eV. Sendo igual ao comprimento total da trajetória do

elétron, alcance R pode ser tomado como um limite superior na determinação de alcances

e, portanto, pode-se concluir que, de acordo com este trabalho, os alcances determinados

por Cole (R

Cole

) e Lea (R

Lea

) e os alcances CSDA1 e CSDA2 para energias menores do que

aproximadamente 1 keV não são boas estimativas para o alcance de elétrons em tecido

uma vez que eles apresentam valores maiores do que aqueles do alcance R.

A razão entre os alcances ponderado e R ( RR

pond

/ ) é igual ao valor médio da fração da

energia do elétron depositada por passo e esta razão (FIGURA 5.4) decresce com o

aumento da energia, isto é, quanto maior a energia do elétron, menor a fração desta energia

depositada em cada passo, o que indica que é no final da trajetória do elétron, onde sua

energia e livre caminho médio são menores, que a maior parte da energia inicial é

depositada.

Considerando a incerteza de ~ 26% associada aos resultados de Davis, os alcances DT, R

e R

se mostraram ser uma boa estimativa para a distância de penetração do elétron em

tecido. Porém, apresentando uma diferença percentual média menor em relação aos valores

de alcance de Davis, o alcance DT é a melhor estimativa, mais especificamente para

elétrons com energias entre 200 eV e 1 keV.

O fato das razões DT / R e R

proj

/ R aumentarem com o aumento da energia inicial do

elétron indica que à medida que a energia aumenta a trajetória do elétron tende a ser menos

tortuosa, ou seja, o elétron tende a apresentar um comportamento mais balístico e menos

difusivo. Porém, o fato destas razões crescerem muito lentamente (aparentando uma

estabilização) a partir de aproximadamente 1 keV indica que, apesar de tender para um

comportamento balístico à medida que a energia aumenta, este dificilmente deixará de ser

em parte difusivo (a não ser, digamos, para energias extremamente altas).

A diferença entre o alcance projetado e o alcance R, para energias entre 200 eV e 50 keV,

é aproximadamente constante e fica em torno de 63%, portanto, a profundidade de

112

penetração de elétrons em tecido, nesta faixa de energia, equivale a aproximadamente 40%

do tamanho da trajetória de elétrons neste meio.

Observando os gráficos das FIGURAS 5.8 e 5.9, nota-se que há uma mudança brusca no

comportamento destas curvas em torno de aproximadamente 1 keV (o que sugere a

ocorrência de uma transição de fase neste ponto). O comportamento dos elétrons com

energias menores do que 1 keV é bastante peculiar e distinto daquele dos elétrons com

energias maiores do que este valor e, portanto, uma maior atenção deve ser dada em

particular a esta importante “classe” de radiação ionizante.

113

APÊNDICE A – Composição do tecido mole (ICRP)

Densidade (g/cm

) = 1,00000E+00

Energia de excitação média I (eV) = 72,300000

Composição

Número

atômico

Peso

(percentual)

1 0,104472

6 0,232190

7 0,024880

8 0,630238

11 0,001130

12 0,000130

15 0,001330

16 0,001990

17 0,001340

19 0,001990

20 0,000230

26 0,000050

30 0,000030

114

APÊNDICE B – Termos, grandezas e unidades relacionados às radiações ionizantes.

1. Grandezas estocásticas e não-estocásticas

Flutuações estatísticas observadas em uma série de medidas de uma determinada grandeza

é um fenômeno bastante familiar à Física. Em muitas destas situações, considera-se apenas

o valor médio do conjunto de medidas e assume-se que, aumentando o número destas

medidas, o valor médio aproxima-se do valor esperado. Este é o caso de grandezas

relacionadas a eventos de deposição de energia e outros eventos discretizados que ocorrem

na interação da radiação com a matéria; efeito que decorre da natureza quântica da própria

matéria e da energia.

Uma grandeza sujeita a flutuações estatísticas é denominada estocástica, e o valor médio

desta, denominada não-estocástica. Assim, a diferença básica entre estas grandezas é que a

grandeza não-estocástica possui um único valor e a grandeza estocástica possui valores que

seguem uma distribuição de probabilidade (ICRU60, 1998).

2. Grandezas e termos relacionados à interação da radiação com a matéria

As grandezas que serão apresentadas a seguir, definidas de acordo com a publicação

ICRU 60 (1998), representam propriedades do meio e do próprio campo de radiação.

Fluência

A fluência ( de partículas), Φ , característica do campo de radiação, é dada por

=Φ , ( B-2.1 )

onde

dN é o número de partículas incidentes em uma esfera com seção transversal de área

da . Esta área deve ser perpendicular à direção de incidência das partículas. A unidade de

Φ é m

-2

Seção de Choque (Cross Section)

A seção de choque,

, de um alvo, para interação com partículas carregadas ou neutras é

, ( B-2.2 )

onde

é a probabilidade de interação com um alvo quando este está sujeito a uma

fluência Φ . A unidade

é m

115

Geometricamente pode-se interpretar

como a área ocupada pelo alvo ou por um disco

em torno deste (ROSSI, 1994): seja um campo de radiação com fluência uniforme e

unidirecional e cujo número de partículas atravessando perpendicularmente um plano de

área unitária é

N . Em uma camada de espessura dl deste material, que contém n alvos

por unidade de volume, o número de alvos por unidade de área é

ndl . Assim, o número de

partículas que sofrem interação,

dN , é:

dlnNdNdN

−

( B-2.3 )

Integrando, tem-se que o número de partículas que não interagiu ao atravessar uma camada

de largura l é:

−

)0(

)(

, ( B-2.4 )

onde )0(N é o número de partículas incidentes. Esta interpretação geométrica de σ não

deve ser aceita literalmente. De maneira geral,

poderia ser definida como a

probabilidade de espalhamento por alvo por unidade de área.

Em estudos de transporte é necessário ainda discriminar os vários processos de

espalhamento através de variáveis adicionais como energia (

h ) e momento ( qh )

transferidos pelo projétil ao alvo no material, ou, equivalentemente, a energia e o ângulo de

espalhamento da partícula espalhada. Formalmente, isto é expresso em termos de seções de

choque diferenciais. Por exemplo,

[

]

hh ddd )( é a probabilidade de ocorrer uma

transferência de energia no intervalo

hd centrado em

h ; e

[

]

)()( )()(

qddqddd hhhh

ωωσ

é a probabilidade de transferência de energia em

hd e

de momento em qd

h centrado em qh . Isto é:

∫

= )(

)()()(

qdd

hhh

( B-2.5 )

∫

= )(

)(

, ( B-2.6 )

onde os limites de integração são determinados por vínculos cinemáticos.

A seção de choque para alvos compostos (moléculas, por exemplo) são usualmente

tratados como se o meio fosse uma mistura de átomos independentes. Isto é válido para a

maioria dos casos, mas ocasionalmente leva a erros, como por exemplo, para fótons de

baixa energia (ROSSI, 1994).

116

Coeficiente de Atenuação

O coeficiente de atenuação,

, para partículas sem carga é dado por

−=

, ( B-2.7 )

onde

N é o número de partículas neutras incidentes e NdN

−

é a fração destas partículas

que interagem ao atravessar uma camada de material de espessura

dl . A unidade de

-1

. O coeficiente de atenuação por massa

, é o coeficiente de atenuação dividido

pela densidade

do material (unidade de m

.kg

-1

), e tem a vantagem de ser independente

da densidade e do estado físico do atenuador. Comparando

com a equação (B-2.3),

tem-se (ROSSI, 1994):

== , ( B-2.8 )

onde

N é o número de Avogadro e

M a massa molar do material.

Coeficiente de transferência de energia por massa

O coeficiente de transferência de energia por massa,

, é dado por

trtr

ρρ

= , ( B-2.9 )

onde

ENdE

é a fração da energia das partículas neutras incidentes que é transferida, sob

a forma de energia cinética, para partículas carregadas do meio de densidade

percorrerem uma distância dl neste.

Coeficiente de absorção de energia por massa

O coeficiente de absorção de energia por massa,

, é

)1( g

tra

−=

, ( B-2.10 )

onde

é a fração da energia das partículas carregadas secundárias que é perdida no

material através de processos radiativos.

Poder de Frenamento (Stopping Power)

O poder de frenamento por massa para partículas carregadas é dado por

117

dES

ρρ

, ( B-2.11 )

onde

dE é a energia perdida por uma partícula carregada ao atravessar uma distância dl

em um material de densidade

. A unidade de

S é J.m

.kg

-1

. S é o poder de

frenamento linear total.

O poder de frenamento por massa total pode ser escrito como a soma de três componentes

independentes:

nucradel

dES

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

ρρρρ

111

, ( B-2.12 )

onde

ρρ

⎟

⎠

⎞

⎜

⎝

⎛

é o poder de frenamento por massa eletrônico (ou de colisão),

devido às colisões com os elétrons atômicos,

rad

ρρ

⎟

⎠

⎞

⎜

⎝

⎛

é o poder de frenamento

por massa radiativo, devido à emissão de radiação de bremsstrahlung nos campos elétricos

do núcleo ou dos elétrons e

nuc

ρρ

⎟

⎠

⎞

⎜

⎝

⎛

o poder de frenamento por massa nuclear,

devido às colisões Coulombianas elásticas em que uma energia de recuo é absorvida pelos

átomos.

Transferência Linear de Energia (Linear Energy Transfer – LET) ou Poder de

Frenamento Linear de Colisão Restrito

A transferência linear de energia,

∆

L , de um material, para partículas carregadas, é

∆

⎟

⎠

⎞

⎜

⎝

⎛

L , ( B-2.13 )

onde

dE é a energia perdida por uma partícula carregada, ao atravessar uma distância dl

no meio, em colisões com elétrons em que a perda de energia é menor do que

∆

. A

unidade de

∆

L é J.m

-1

∞

L é igual a

col

S , o poder de frenamento linear total devido a

colisões.

118

APÊNDICE C – Cinemática de colisão

1. Cinemática do processo de espalhamento

A cinemática lida com as conseqüências da conservação de energia e momento nos

processos de colisão. A FIGURA C-1.1 ilustra esquematicamente uma colisão de dois

corpos e introduz algumas notações para caracterização de cada partícula:

= energia cinética

= momento

v = velocidade

c = velocidade da luz no vácuo

m = massa de repouso

mcTE += = energia total

cv=

212

)1(

−

−=

βγ

Algumas dessas grandezas se relacionam como a seguir:

22422

)1(

cpcmmc

E −==

−

( C-1.1 )

βγ

TmcT

cmE

p =

−

22422

( C-1.2 )

Quando

mcT << pode-se usar a cinemática não-relativística.

Muitos cálculos de cinemática são simplificados se a colisão é tratada no sistema de

referência do centro de massa (CM); ele é definido de forma que o momento total das

partículas – antes e depois da colisão – é nulo (FIGURA C-1.1). O sistema de referência

em repouso em relação ao laboratório recebe este nome (L); tipicamente, mas não sempre,

em L o alvo está em repouso.

Classicamente, a velocidade do CM em relação a L é (FIGURA C-1.1)

= . ( C-1.3 )

Aqui, os índices 1 e 2 referem-se ao projétil e ao alvo respectivamente. As grandezas em L

são escritas em letras minúsculas e, em CM, em maiúsculas. As velocidades em L e CM

podem ser obtidas uma a partir da outra:

119

= ( C-1.4 )

−= . ( C-1.5 )

Na colisão de um projétil com um alvo muito pesado (

< ) não há distinção entre L e

CM uma vez que

v é praticamente nula. Este é o caso, por exemplo, de um elétron

espalhado por um átomo (

−

≅mm ).

Em colisões elásticas

4324321

; TTTT mmmm

+ . ( C-1.6 )

Em colisões inelásticas a energia interna das partículas é envolvida nos processos de

transferência de energia, isto é, ocorrem conversões de energia entre massa e energia

cinética. Daí segue que

QTTTT Qmmmm

−

4324321

; ( C-1.7 )

e a reação é chamada exotérmica (endotérmica) se 0>Q (0

Q ).

O número de variáveis ( , , ,

etc) que podem ser obtidas a partir da cinemática é igual

ao número de equações que expressam, independentemente, conservação de energia e

momento. Para a colisão representada na FIGURA C-1.1 são três as equações:

QTTT

−

( C-1.8 )

44331

coscos

ppp

( C-1.9 )

4433

sensen0

. ( C-1.10 )

A energia transferida na colisão é

FIGURA C-1.1: Espalhamento de dois corpos no referencial do laboratório L

(esquerda) e no referencial do centro de massa CM (direita).

120

h , ( C-1.11 )

onde considerou-se

mm = e

(mesmas partículas antes e depois do

espalhamento).

O valor máximo de

h ocorre quando 0

e 180

(colisão frontal):

max

)(

h . ( C-1.12 )

Quando

mm << ,

()

121max

4 Tmm=

h , o que significa que a transferência de energia em

colisões elásticas é extremamente pequena. Este é o caso do espalhamento elétron-átomo,

elétron-molécula ou elétron-sólido.

Finalmente, em reações em que 0

Q , é necessário que o projétil tenha uma certa energia

mínima para que a reação ocorra. Esta energia é dada por:

[

]

min

)(2

cmcmQQ

++−−

= ( C-1.13 )

(ROSSI, 1994).

121

REFERÊNCIAS BIBLIOGRÁFICAS

ATTIX, F.H.; ROESCH, W.C.; TOCHILIN, E. Radiation dosimetry: fundamentals. 2.ed.

New York, Academic Press, 1968, v.1.

BERGER, M.J. Monte Carlo calculation of the penetration and diffusion of fast charged

particles, in

Methods in Computational Physics, vol. 1, p. 135-215, eds. ALDER, B.,

FERNBACH, S.; ROTENBERG, M., Academic Press, New York, 1963.

BERGER, M.J.; COURSEY, J.S.; ZUCKER, M.A. ESTAR, PSTAR, and ASTAR:

Computer Programs for Calculating Stopping-Power and Range Tables for Electrons,

Protons, and Helium Ions. Version 1.2.2. Disponível em: http://physics.nist.gov/Star

Acesso em 15.05.05. National Institute of Standards and Technology, Gaithersburg, MD,

2000.

BETHE, H.A. Zur Theorie des Durchgangs scheneller Korpurkularstrahlen durch Materie.

Ann. Physik, vol. 5, p. 325-400, 1930.

BETHE, H.A. Bremsformel für Elektronen relativistischer Geschwindigkeit.

Z. Physik,

vol. 76, p. 293-299, 1932.

BIELAJEW, A.F.; ROGERS, D.W. Variance-reduction techniques, in

Monte Carlo

Transport of Electrons and Photons, eds. JENKINS, T.M.; NELSON, W.R.; RINDI, A.

Plenum, New York, p. 407-419, 1988.

BLOOMER, W.D.; ADELSTEIN, S.J.

Pathobiological Annuary. 8, 1978.

CHARLTON, D.E.; BOOZ, J.; FIDORRA, J.; SMIT, TH.; FEINENDEGEN, L. In:

PROCEEDINGS OF SIXTY SYMPOSIUM OF MICRODOSIMETRY, Harvard

Academic Publ. Ltd., London, 1978.

COLE, A. Absorption of 20 eV to 50.000 eV Electron Beams in Air and Plastic.

Radiation

Research, vol. 38, p. 7-33, 1969.

DAVIS, M. Range Measurement of Low-Voltage Electrons.

Phys. Rev., vol. 94, n. 2,

p. 243-245, 1954.

EMFIETZOGLOU, D.; KARAVA, K.; PAPAMICHAEL, G.; MOSCOVITCH, M. Monte

Carlo simulation of the energy loss of low-energy electrons in liquid water.

Physics in

Medicine and Biology. vol. 48, p. 2355-2371, 2003.

122

FANO, U. Penetration of protons, alpha particles and mesons.

Ann. Ver. Nucl. Sci.,

vol. 13, p. 1-66, 1963.

FERMI, E. Nuclear Physics University of Chicago Press, 1950.

FERNÁNDEZ-VAREA, J.M.; MAYOL, R.; BARÓ, J.; SALVAT, F. On the theory and

simulation of multiple elastic scattering of electrons.

Nucl. Instrum. Meth. B, vol. 73,

p. 447-473, 1993.

GOUDSMIT. S.; SAUNDERSON, J.L. Multiple scattering of electrons.

Phys. Ver.

vol. 57, p. 1459-1464, 1940.

HUBBELL, J.G. Mass attenuation and energy-absorbing coefficients.

Int. Appl. Radiat.

Ist. vol. 33, p. 11, 1982.

HUMM, J.L.; ROESKE, J.C.; FISHER, D.R.; CHEN, G.T.Y. Microdosimetric concepts in

radioimmunotherapy

. Medical Physics, vol. 20, p. 535-541, 1993

INTERNATIONAL COMMISSION ON RADIATION UNITS AND MEASUREMENTS

(ICRU).

Radiation Quantities and Units, ICRU Rep. Nº. 33, ICRU, Bethesda, Maryland,

1980, (ICRU, 33).

INTERNATIONAL COMMISSION ON RADIATION UNITS AND MEASUREMENTS

(ICRU).

Stopping Powers for Electrons and Positrons, ICRU Rep. Nº. 37, ICRU,

Bethesda, Maryland, 1984, (ICRU, 37).

INTERNATIONAL COMMISSION ON RADIATION UNITS AND MEASUREMENTS

(ICRU).

Fundamental Quantities and Units for Ionizing Radiation, ICRU Rep. Nº. 60,

ICRU, Bethesda, Maryland, 1998, (ICRU, 60).

JAMES, F. Monte Carlo theory and pratice.

Rep. Prog. Phys., vol. 43, p. 1145-1189, 1980.

KASSIS, A.I.; ADELSTEIN, S.J.; HAYDOCK, C. et al.

Radiation Research, 84, 1980.

LEA, D.E.

The Action of Radiation on Living Cells. Cambridge University Press,

Cambridge, p. 23, 1947.

LÉCUYER, P. Efficient and portable combined random number generators.

Commun.

ACM, vol. 31, p. 742, 1988.

MALAMUT, C.

Difusão de elétrons em meios simuladores de tecidos. 81p. Tese

(Doutorado em Física). Pontifícia Universidade Católica do Rio de Janeiro, Departamento

de Física, 1987.

123

MALAMUT, C.; PAES-LEME, P.J.; PASCHOA, A.S. Diffusion of low-energy electrons

in tissue-like liquids.

Radiation Research, vol. 132, p. 168-177, 1992.

MAUSNER, L.F.; SRIVASTAVA, Selection of radionuclides for radioimmunotherapy.

Medical Physics, vol. 20, n. 2, p. 503-509, 1993.

MEESUNGNOEN, J.; JAY-GERIN, J-P.; FILALI-MOUHIM, A.; MANKHETKORN, S.

Low-energy electron penetration range in liquid water.

Radiation Research, vol. 158,

p. 657-660, 2002.

MOLIÈRE, G.Z. Theorie der Streuung scheneller geladener Teilchen II Nehrfach und

Vielfachstreuung.

Z. Naturforsch, A 3A, p. 78-97, 1948.

PENELOPE - A CODE SYSTEM FOR MONTE CARLO SIMULATION OF

ELECTRON AND PHOTON TRANSPORT. AEN, NEA, 2001.

PRATT, R.H.; TSENG, H.K.; LEE, C.M.; KISSEL, L.; MacCALLUN, C.; RILEY, M.

Bremsstrahlung energy spectra from electrons of kinetic energy 1 keV < T

< 2000 keV

incident on neutral atoms 2 < Z < 92,

Atomic Data Nucl. Data Tables, vol. 20, p. 175,

1977. Errata in vol. 26, p. 477, 1981.

ROGERS, D.W.O; BIELAJEW, A.F. Monte Carlo Techniques of Electron and Photon

Transport for Radiation Dosimetry, in

The Dosimetry of Ionizing Radiation. vol. 3, cap. 5,

p. 427-539, eds. ATTIX, F.H.; BJARNGARD, B.E.; KASE, K.R. Academic Press, 1990.

ROSSI, B.

High-energy particles. New York: Prentice-Hall, Inc., 1952.

ROSSI, H.H; ZAIDER, M. Microdosimetry and its applications. New York: Springer,

1994.

RUBINSTEIN, R.Y.

Simulation and the Monte Carlo Method. Wiley, New York, 1981.

SALVAT, F. Algorithms for random sampling from single-variate distributions.

Comput.

Phys. Commun., vol. 46, p. 427-436, 1987.

SALVAT, F.; FERNANDEZ-VAREA, J.M.; BARO, J; SEMPAU, J.

An Algorithm and

Computer Code for Monte Carlo Simulation of Electron-Photon Showers. Informes

Tecnicos CIEMAT no 799, Madrid: CIEMAT, 1996.

SALVAT, F.; FERNÁNDEZ-VAREA, J.M.; ACOSTA, E.; SEMPAU, J. Penelope: A

Code System for Monte Carlo Simulation of Electron and Photon Transport. In:

WORKSHOP PROCEEDINGS, 5-7 Nov. 2001, Issy-Les-Moulineaux.

Proceedings... Issy-

Les-Moulineaux: Nuclear Energy Agency, 2001. 234p.

124

SALVAT, F.; FERNÁNDEZ-VAREA, J.M.; ACOSTA, E.; SEMPAU, J. Penelope:

Manual.txt. Version 2003. Disponível em:

http://geology.wisc.edu/~johnf/g777/777Penelope/manual.pdf

. Acesso em 11.12.03.

SOBOL, I. O Método de Monte Carlo

. Moscou: Mir, 1983.

STERNHEIMER, R.M. The density effect for the ionization loss in various materials.

Phys. Rev., vol. 88, p. 851, 1952.

STERNHEIMER, R.M.; SELTZER, S.M.; BERGER, M.J. Density effect for the ionization

loss of charged particles.

Phys. Rev. B, vol. 26, p. 6067, 1982.

STEWART, R.D.; WILSON, W.E.; McDONALD, J.C.; STROM, D.J. Microdosimetric

properties of ionizing electrons in water: a teste of the PENELOPE code system.

Physics

in Medicine and Biology, vol. 47, p. 79-88, 2002.

TAUHATA, L.; SALATI, I.P.A.; DI PRINZIO, R.; DI PRINZIO, A.R.

Radioproteção e

dosimentria: fundamentos. Rio de Janeiro: Instituto de Radioproteção e Dosimetria, 2002.

TISLJAR-LENTRILIS, G.; HENNEBERG, P.; MIELKE, TH.; FEINENDEGEN, L.E. In:

PROCEEDINGS OF THE SIXTY SYMPOSIUM ON MICRODOSIMETRY, Brussels,

Comission of the European Communities, Harwood, London, 1978.

TURNER, J.E.

Atoms, Radiation, and Radiation Protection. 2.ed. New York, John Wiley

& Sons, Inc., 1995.

WILLIAMS, E.J.

Proc. Roy. Soc. Londres, A130, p. 310, 1930

FIGURA 1.1:

Absorção de

elétrons

monoenergétic

os.

R é

Livros Grátis
( http://www.livrosgratis.com.br )
 
Milhares de Livros para Download:
 
Baixar livros de Administração
Baixar livros de Agronomia
Baixar livros de Arquitetura
Baixar livros de Artes
Baixar livros de Astronomia
Baixar livros de Biologia Geral
Baixar livros de Ciência da Computação
Baixar livros de Ciência da Informação
Baixar livros de Ciência Política
Baixar livros de Ciências da Saúde
Baixar livros de Comunicação
Baixar livros do Conselho Nacional de Educação - CNE
Baixar livros de Defesa civil
Baixar livros de Direito
Baixar livros de Direitos humanos
Baixar livros de Economia
Baixar livros de Economia Doméstica
Baixar livros de Educação
Baixar livros de Educação - Trânsito
Baixar livros de Educação Física
Baixar livros de Engenharia Aeroespacial
Baixar livros de Farmácia
Baixar livros de Filosofia
Baixar livros de Física
Baixar livros de Geociências
Baixar livros de Geografia
Baixar livros de História
Baixar livros de Línguas

Baixar livros de Literatura
Baixar livros de Literatura de Cordel
Baixar livros de Literatura Infantil
Baixar livros de Matemática
Baixar livros de Medicina
Baixar livros de Medicina Veterinária
Baixar livros de Meio Ambiente
Baixar livros de Meteorologia
Baixar Monografias e TCC
Baixar livros Multidisciplinar
Baixar livros de Música
Baixar livros de Psicologia
Baixar livros de Química
Baixar livros de Saúde Coletiva
Baixar livros de Serviço Social
Baixar livros de Sociologia
Baixar livros de Teologia
Baixar livros de Trabalho
Baixar livros de Turismo
 
 