Barichello, Liliane Basso

Reichert, Janice Teresinha

Universidade Federal do Rio Grande do Sul

Solução em forma fechada para a equação íntegro diferencial de transporte de fótons

Lume 5.8

155N 1984-820X ... -.. ~ ffi @rn[;ffi@ • ... ····~ e CNMAC2010 • :. • ••• XXXII I CONGRESSO NACIONAL •• ....... DE MATEMÁTICA Af>UCADA \!'Jb® PUI!ILICACÁO DA Sl!l MAC ••••• E COMPUTACIONAL Solução em forma fechada para a equação íntegro diferencial de transporte de fótons Janice Teresinha Reichert, Departarrento de Matemática, Universidade Tecno1ógica Federal do Paraná 85504-300, Pato Branco, PR, E-mail: j aruce@utfPr.edu. br, Liliane Basso Barichello Instituto de Matemitica, U niversldade Federal do Rio Grande do Sul 91509-900, Porto A1egre, RS, E-mail: lbaric@mat. ufrgs. br Resumo: Neste trabalho, uma solução em forma fechada é proposta para a equação íntegro diferencial de transporte de fótons unidimensional. A dependência energética no problema é abordada de forma discretizada. O problema resultante, com dependência angular e espacial, é resolvido pelo método analítico de ordenadas discretas (ADO), resultando em uma solução analítica em termos da variável espacial. Simulações numéricas avaliam o chamado fator de buildup, usado como fator de correção em procedimentos práticos associados ao cálculo de doses em tratamentos de radioterapia. Palavras-chave: Equação de transporte, núcleo de Klein-Nishina; radioterapia; método analítico de ordenadas discretas. 1 Introdução Uma das quantidades de interesse em tratamentos por radioterapia é a energia depositada, ou dose de radiação. Para corrpreender corro ocorre a interação da radiação com a matéria, e analisar a energia depositada, fuz-se necessário um bom entendirrento da distribuição do fluxo (intens iiade) de radiação num determinado rreb [5], que pode ser obtido através da soh.!ção da equação íntegro diferencial de transporte de partícu1as, parti:u1armente de fõtons [ 6]. Para o caso de baixas energias, a caracteristi:a do chamado espafuarrento Corq:>ton, predominante no processo [6, 9], introduz uma dificuldade maior no tratarrento do terrrn integral da equação. Resuhados clássicos são conhecidos [9] para problemas unidirrensiona:is em rreio infinito, utilizando rrétodo conhecido corro o rretodo dos morrentos, que se resume na sohlção de uma equação integral de Voherra por aproximação n~rica. A sohlção do prob-Jema emrreb finito, onde esta abordagem citada anteriorrrente não se aplica, é serrpre um desafio e, em georretr:ias corrplexas é em geral tratato por abordagens probabilisti:as [4]. Uma dificu1dade adicbnal do prob1ema, além da conheciia corq:>Jexidade dos rrodeb s que envolvem a equação de Bohzmann linear, se encontra na dependência dos limites de integração da equação à var:iável corq:>rimento de onda. Neste trabafuo, na seção 2, introduzirms a equação de transporte para fótons com núc1eo de es-pa1hamento de K1ein-Nishina, que define o prob1ema básico a ser considerado. A seguir, na seção 3, proporros um tratarrento discret:izado para a variável energia e apresentamos o rrétodo ADO [ 1] para tratarrento da dependência angu]ar da equação, resuhando em uma sohlção em forma fechada para o prob1ema, analítica em terrms da variável espacial FinaJrrente, nas seções 4 e 5, apresentarros resuha-dos nurréri:os e considerações finais . 493 ISSN 1984-820X ... -.. ~ffi@]~ • ..... ··~ e CNMAC2010 • :. • ••• XXXII I CONGRESSO NACIONAL ...... ;.. DE MATEMÁTICA APLICAOA li_'go® PUISLICACAO DA SISMAC ••••• E COMPUTACIONAL 2 Equação de Transporte para Fótons Cons:i:ieremos a equação de transporte para iõtons escrita na furma [6], ô ~ 2[- 1 11 -8 J(.r.À.. fl ) ,- a (>.. )J(.r. À. fl) =a --1~ (/L )x f 2 1=0 H (>..' . >..)JJ1(1 + X - >.. )1), (!/ )J(.r . >..' . f/ )df/d>..' . (1) max{:lo . .X- 2} - J onde 1 (.r. À. 11 ) é a n tensidade de radiação, ~~: E (O .. r0 ) representa a variável espacial , 1 ~ E r - 1. 1] é o cosseno do angub polar, À E r>..o, ) é o comprimento de onda e Ào O comprimento de onda associado à energia do :teix.e n c idente de radiação, no domínio. A lém disto, O'(À) representa o coeficiente. de atenuação fuear. Considera-se, o núcleo de K1en-N ishina para espahamento Compton [2, 6, 7, 9], A'(>..'.>.. ), presente na equação (1), defin:i:io como sendo nulo, exceto se>..' :::; >.. :::; X + 2, sendo, neste caso, dado por (2) A:irxia, = l\'A ZI'YTP =O 10061 Zp r · 1] n A . A. c llt . (3) onde A é a massa atônrica, Z é o número atômico, .\',-~ é o número de Avogadro, o-r é a seção de choque de Thomson, e p é a densidade do meio. Devido à característica do espalliamento Compton, e confOrme defini;ão do nú:leo de K Jein-Nishina, o tenro iltegral na equação (1), de :làto, tem limite interior igual a >..0 ou >.. - 2, o que fur maior, confOrme já fui observado por [2]. Neste trabafuo, propõe-se urna solução de caráter analíti::o, para a aprox:inação em ordenadas discretas da equação (1), considerando que a n tensidade de radiação 1 (.r. À.. fl ), seja con-hec:i:ia nas fronteiras do dominio: placa plana de espessura .r0 . 3 Desenvolvimento Cons:i:ierando que o termo integral na equação (1) tem li:níes dependentes da variável comprimento de onda >.., o que dificulta a proposi;ão de um esquema de quadratura, propomos neste trabaho, um esquema que assegure que o valor nw.:.c{ >.. 0 , À - 2} seja um dos nós da quadratura, caso contráril o problema discret:izado levar :R a um sistema com número maior de incógní as do que de equações. As-sim, primeiramente, fixa-se um valor máximo, denotado por /\ 1.· , como aproximação numérica do limite superior (infirrito) do intervalo de definição da variivel de n teresse, À E r>..o. ). Posterilrmente, o espaçamento unifOrme entre os pontos, h= (>..F- /\ 0 )/n é definido baseado na escofua de um número de :interva1os 11 apropriado para que a propriedade acima mencionada se cumpra. A:irxia, >..o = 0.511/ Eo é um valor conhecido. Assin, J l = n- 1 será o número de intervalos (grupos) de energia. Fazendo >.. 1 = >..0 o prmeiro valor discret:izado, os demais valores para a discret:ização da variivel energ:B são cakuJados como, À; = >..;_ 1 - h, i= 2. 3 ..... Jl, garantindo que ma.r{/\1 . À;- 2}, i = 1. 2 .... Jl esteja sempre entre os valores discretizados de >.. . Para cada À;, i = 1. .... J l , resolve-se então a equação de transporte, ô ~ 2[- 1 11-8 l (.r. À.; . fi ),- O' (À;) l (.r. À; . fL ) =a --JH11) f 2 1=0 X }{ (>..' . À; )F1(1 + À.1 - À; )JJ1(!/)J (.r . À 1 . f /)df/ d>..'. ( 4) ma.,-{ .X 1 . .Xi - 2} - J A :integral em >..' será resolvidapelt regra dos trapézios . Cl::lamanio, 494 ISSN 1984-820X ... -.. ~ cili; @IT][IJ2@ • ..... ··~ e CNMAC2010 • :. • ••• XXXIII CONGRESSO NACIONAl ...... ; •• DE MATEMÁTICA APLICADA ~o® PUISLICAÇAO OA SISMAC ••••• E COMPUTACIONAl (5) tem-se, ü ~ 2[-l fi-nl (.r . À1. fi )- a (>..;) l (.r. À.; . fl ) =a --J)I(ft) Uf 2 1=0 X ]( (>..' . À;)l)l( l + X - À, )l H f/ )J(.r. >..' . f/ )df/d>..' . (6) .>.,. - 1 e usando a regra dos trapézi>s, 1~ ~,L (.r. À;. fi ) + a (>..; ).i ( r, >..u~) =a~ 2l; 1 F1 (t~)(À.; ; )..,. )[1\. ( >.. ,., À; )F1(1 + >.., - >..J x l = O 11(f/ )1(.r . À1 . ft' )df11 + }( (>.. ;. À. ;)FI (l + À;- À;) 11(f/ ).i (. r . À. ;. f t')dft' ]. (7) I onde >.., é um valor fuco (para cada i) . Chamando 1(x . À; . ft ) = 1;(.r. ft ) . e denotan:io por C 1, a matriz de ordem J1 x J1 com coefi:ientes não rruhs definidos por, (c) - (>..;- >..,) 1· 1) ( I ' ' ) I i .r -2 \ ri I + I\ r - 1\j · pode-se escrever a equação (7), u tilizando a notação matricial, D L n ~ fi -D I (.1·, f ~) + L I (.1:, !L) = - (2l + l )11{t~ )C 1 .r 2 1= 0 11(fi.' )l (. r, fi.' )dt/, I (8) (9) (10) (11) para .r E (O .. ro ) e fi. E: [- 1. 1], onde L é a matriz diagonal com elementos a; = a (>..;), i = 1, 2 .... i\1, e, I (.r. f~)= (JI(:r . f ~) . h (.t . f ~) ..... .I,,,(.r.fi ))"' . (12) define o vetor J\1 x 1 com componentes 1 (.r. À. ; . f~ ) . No prob ~ma resu1tante o tratamento da variável angular é reito pelo ~todo ana 1ítico de ordenadas discretas (ADO) [1]. 3.1 Tratamento da variável angular pelo método de ordenadas dhcretas A equação (11) em pode ser escrita em termos de ordenadas discretas, d (\ ~ ~ fi ; d.r I (.r.ft;) - L I (.:r. f1;) = ? (21 + l )J>1 (fi ; )C 1 w., I1,.Jr ). - 1=0 . - 1 e L N (13) d n ~ 1 ~ - ti;- l (.r. - ti;) + L I (.1:. -f1;)=- (-1)(2l + l).IJI (I' i)C I u•., I I J r ). (14) d.r 2 ' 1= 0 ,,= I 495 ISSN 1984-820X ... -.. ~ cili; @IT][IJ2@ • ..... ··~ e CNMAC2010 • :. • ••• XXXIII CONGRESSO NACIONAl ...... ; •• DE MATEMÁTICA APLICADA ~o® PUISLICAÇAO OA SISMAC ••••• E COMPUTACIONAl para i = l. · · · . . \ ' , e onde (15) Aqu~ denotarros por { Jt., 1 os pontos de quadratura e { 11\1 os pesos da quadratura definidos no intervalo '0 . 1' . Buscamo soh.Jções da fOrma exponencial para as equações (13) e (14), escrevemos (16) substíuindo esta expressão nas Eqs. (13) aro (14) e após algumas manipulações algébri:as [8], podemos escrever, para i = 1. · · · . .\1, o sistema on:ie <l>t .. , = 11(Jt.,) [<t> (v. Jt., ) + ( - 1)1<t> {v. - Jt., )]. e f é a matriz identnade .\1 x .\1. Introduzindo os vetores de ordem .\1. \', e r T T 'f' ' '/' 4> (1/)= 4> (v. -JLJ) . <t> (1/. -Jt2) ..... 4> (1/. - Jt.v ) . on:ie T denota a operação transposta, as matrizes quadradas (.\1.\' x .\1.\' ), W = dia.g{ .... u,,. ... f . M = dío.g{ .... ,t,l .... } . D = dío.g { .... L , ... } . e as matrizes de ordem ( J\1 N x J\1) as Eqs. (17) e (18) podem ser escd:as como, e L 1 n ~ 1 (D ~ - M )<t>_(v) = - (-1) (2t -l)n ,c1G1(1/) . v 2 1=0 com Agora, chamando (19) (20) (21) (22) (23) (24) (25) (26) (27) (28) (29) (30) 496i ISSN 1984-820X ... -.. ~ cili; @IT][IJ2@ • ..... ··~ e CNMAC2010 • :. • ••• XXXIII CONGRESSO NACIONAl ...... ; •• DE MATEMÁTICA APLICADA ~o® PUISLICAÇAO OA SISMAC ••••• E COMPUTACIONAl fàzendo a soma e a d iferença das Eqs. (26) e (27) obtemos EX = ~ Y e H Y = ~X . v v o me L E = (D - ~ ~ (2l- 1 )n ,C / 1 + (-1)1'n rw )M I 1=0 L n ~ 11 T I H =(D - ? (2l + 1)n ,C,[1 -(- 1) n 1 W )M - 1= 0 X = M U . Y = MV . Utilizando as equações (31) obtém-se do is problemas de autovalores, orxle .-\ = 1/ 12. (H E)X =/\X (EH )Y = /\Y (31) (32) (33) (34) (35) (36) Os autova hres /\J e autovetores correspondentes X (.-\ J ) , j = 1. 2 . .. .. .\1 .\ ' são utilizados para esc r-ever 1 I A <l> +(vj) = 2M - ( I - v1E )X (Ãj) 1 I A <l> _ (vj) = 2M - (I - t11E )X (ÃJ) · para j = 1. 2 . .. . ( .\1 .,.. 1 )X A expressão f é a matriz :identidade de ordem .\1 .\' x .\1 X Agora definindo T 1' T T l + (.r )= [l (X. flJ) . I (.r. f l2 ) ..... 1 (.t. ft .v)] . T T '/' '/' (.r)= [l {:1:, - fl J). I (.r . - tl2) . ... . 1 (.r . - fi N) ] . (37) (38) (39) (40) podemos escrever, em forma :fechada, a solução em ordenadas discretas, para o prob e ma homogêneo dado pelas Eqs. (13) e (14) como AI N I + (.r ) = ~ [Aj <l> +(vj)e- :c/111 + llj <l> (vj)<' (.ru .r)/ ".1 ]. (41) J = l e .\IX (.r)= ~ [Aj <l> _(v1)e-4"J - HJ<l> +(vj)c (•·o •·l/".1 ] (42) j=l orxle as constantes AJ e BJ são detenninadas impondo-se as cond i';ões de contorno. 4 Resultados numéricos e aspectos computacionais O cák u1o da dose de radiação, uma das quantidades de ilteresse que podem ser avaliadas a partir da so]ução desenvolvida nas seções anterDres, pode ser :feita, segurxlo [2, 7], como (43) 497 155N 1984-820X ... -.. ~ ffi @rn[;ffi@ • ... ····~ e CNMAC2010 • :. • ••• XXXII I CONGRESSO NACIONAL •• ....... DE MATEMÁTICA Af>UCADA \!'Jb® PUI!ILICACÁO DA Sl!l MAC ••••• E COMPUTACIONAL onde o fluxo <Í? (p.0:c, À) é definiio por, ,,1 <Í?(;1.o:c, À) = I (Jl.o::c, À, ,,,)cl;,,, (44) 1 mcc2 = 0.511/ E0, sendo E0 igual a energia incidente, Jl.o o coeficiente de atenuação do rreio, para a energia incidente E0 e ' ''a(À) o coeficiente de absorção, utilizado aqui corro sendo o do ar. Para avaliação n~ca da equação (43), o rresrm esquema de discretização e quadratura descrito na seção (3) é utilizado. O fator de buildup, que é norrnahrente usado corro tnn parâmetro de correção (para evitar a sohlção da equação que envolve o termo integral de espalharrento) entre o fluxo que inchri o espafuarrento e o fluxo deviio aos fótons não espafuados [3 ], é caku1ado corro B( ") _ Dy(JI.o~: ) Jl.Q~L - ( ) , Do Jl.o ~c (45) onde a expressão no denominador é a dose devido à radiação primllia, ou seja, ao fluxo de radiação não colidido <Í?( Jl.o:~;, À o) . Para o problema teste aqui abordado definirms, corro condições de contorno, I (O, Ào, ,,,) = 1, para ,,, > O, I (O, Àí, J',) = O, para ,,, > O, i = 1, 2, ... , M, I(:~;o , Àí, ,,,) = O para ,,, < O, i = 1, 2, ... , M , (46) (47) (48) e consiierarms uma placa p1ana horrogênea de água (espessura definida emterrms de ,,.0:~; ), com energia incidente de 1M e v, sendo, neste caso, Jl.o = 0.070G o va1or do coeficiente de atenuação. Usamos nove grupos de energia para a d:i'5cret:ização da variável COIIJ>rirrento de onda e N = 20 para o rrétodo ADO. Detenninamos o fator de builiup, e coiiJ>ararms com resuhados obtidos com outra abordagem, que considera a dependência espectral corro contínua, apresentada em [8], Tabe1a 1: CoiiJ>aração do futor de builiup. ,,.0 :~; Este trabafuo Ref [8] 1 2.5281 2.1111 2 3 4 5 6 2.8933 3.2588 3.6569 4. 1116 4.6399 2.8704 3.5124 4.0692 4.5481 4.9513 Observarms uma concordância minimarrente satisfatória entre as duas abordagens utilizadas. Al-gurms diferenças, do caso aqui apresentado e na abordagem contínua, podem ser atribuídas a resohlção tot:ahrente analítica utilizada no caso contínuo (baseada em uma expansão espectral na variável com pr:iJ:rento de onda) corrparada com a abordagem aqui utilizada através da regra dos trapézi:>s para a resohlção das integrais. No entanto, a proposta aqui descrita pode ser considerada mais s:irrples e pos-sivehrente mais fucihnente inserida no tratarrento de problemas ma:i'5 coiiJ>lexos. De acordo com o nosso conhecirrento não existem na literatura outros resuhados n~cos para este prob1ema aplicado a p hcas planas. Na tabeh seguinte apresentarms uma descrição de va1ores para diferentes ordens de aproximação, para esta abordagem, considerando o meio corro água, com energia incidente de 1M ev e nove grupos de energia, em :fimção do núrrero de pontos de quadratura considerados para o rrttodo ADO, 498 155N 1984-820X ... -.. ~ ffi @rn[;ffi@ • ... ····~ e CNMAC2010 • :. • ••• XXXII I CONGRESSO NACIONAL •• ....... DE MATEMÁTICA Af>UCADA \!'Jb® PUI!ILICACÁO DA Sl!l MAC ••••• E COMPUTACIONAL Tabela 2 : Convergência numéri:a do método ADO. N tJ-o:r: = 1 6 2.5291 10 2.5283 14 2.5282 20 2.5281 24 2.5281 28 2.5281 Como observa-se na tabela 2, com N = 20 obtém-se a concordância de quatro dígitos significativos, o que indica uma ráp ida cotNergênc it numérica do método. Os resultados numéricos são obtidos com baixo es fOrço computacional e gerados utilizando a lin-guagem FOR1RAN e executadas em um notebook Centrino 2 Duo 2GHz, em menos de 5 segundos. O bom desempenho computacionaL deve-se principalmente à utilização do método analítico de orde-nadas discretas (ADO) que reduz o prob lema original ao cá kulo de um esquema de quadratura no seml- interva b , que implica num ganho de precisão na discretização em relação à esquemas clássicos de quadratura utilizados em problemas de transporte. 5 Comentários Finais Apresentou-se uma solução em furrna fechada, analítica em termos da variáve l espac itL para a equação de transporte de fõtons com espalhamento Cornpton Foi usada uma abordagem d iscreta da dependênc it espectraL para obtenção de resuhados numéri:os, para o problema da placa plana submetida a condições de contorno de radiação incidente conhecidas nas extremilades . A abordagem de quadratura proposta permitiu um tratamento satisfàtório ao problema de dependência dos limites de integração em re ltção à variável comprimento de onda. Referências [1] L.B. Barichelb , C. E. S iewert, A Discrete-Ordinates Solut :ion fura Non-Grey Model With Com plete Frequency Red istribution, Quantitative Spectroscopy & Radi.ative Transfer ,62 (1999), 665-675. [2] U. Fano, L. V. Spencer, M.J. Berger, ' 'Encycloped it of Physics" - Vo h.Jme XXXVIII/2, Neutrons and Re ltted Gamma Ray Problems, Springer-Verlag, Berlin, 1959. [3] H. Goldstein, ' 'Fundamental Aspects ofReactor Shielding", Addison-Wesley Publishing Company, INC.,1959. [4] H. Hirayama, K. Shin, App lication of the EGS4 Monte Carlo Code to a Study of Muhilayer Gamma-Ray Exposure Builiup Factor of up to 40 mfP, Journal ofNuclear Science and Technology, 35 (1998), 816-829. [5] H.E. Johns, J.R. Cunninghan, ''The Physi:s ofRad :iology", Thomas Springfie ll Press, 1983. [6] J. K. Shultis, R. E. Faw, ' 'Radiation Shielding", American Nuc 1ear Society, Inc., 2000, La Grange Park, Illino is, USA. [7] J. J. Duderstadt, W. R. Martin, ' Transport Theory", John Wiley & Sons, Inc., 1979. [8] J.T. Reichert, '~bordagens Ana td:i:as para Problemas de Transporte de Radiação com Dependênc it Espectraf ', Tese de Doutorado, PROMEC-UFRGS, 2009. [9] J. Wood, "Cornputational Methods in Reactor Shielding", Oxford, Pergamon Press,1982. 499 

Portuguese

https://lume.ufrgs.br/bitstream/10183/181325/1/000772876.pdf

Solução em forma fechada para a equação íntegro diferencial de transporte de fótons

Abstract

Similar works

Full text

Available Versions

Universidade Federal do Rio Grande do Sul

Lume 5.8