The application of finite element method in space dependent neutron diffusion equation is shown in this paper The diffusion equation is first discretized using rectangular elements Applying piecewise polynomials in a variational form of the diffusion equation a system of linear algebraic equations is obtained Choleski decomposition is used to factorise the global matrix and numeration of nodes and finite elements is performed by line dissection The multigroup two-dimensional diffusion code KONEL was written based on this method (author)Society for Electronics,Telecommunications, Automation, and Nuclear Engineering; Belgrade (Yugoslavia); ETAN '86: 30. Conference; 30. Konferencija za ETAN; Herceg Novi (Yugoslavia); 2-6 Jun 1986

Petrović, Ivan M.

Croatian

VinaR - Repository of the Vinča Institute of Nuclear Sciences

IX.16 9 
XXX JUGOSLOVENSKA KONFERENCIJA ETAN-a. HERCEG-NOVI, 2—6. . JUNA 1986. 
Ivan Petrović 
Institut za nuklearne nauke "Boris Kidrič" - Vinča 
OOUR Institut za nuklearnu energetiku i tehničku fiziku "NET" 
PRIMENA METODE KONAČNIH ELEMENATA NA REŠAVANJE VIŠEGRUPNE 
DVODIMENZIONE DIFUZIONE JEDNAČINE X RAČUNARSKI PROGRAM KONEL 
THE APPLICATION OF THE FINITE ELEMENT METHOD TO THE MULTIGROUP 
TWO-DIMENSIONAL NEUTRON DIFFUSION EQUATION AND THE KONEL CODE 
SADRŽAJ - Analizirana je primena metode konačnih elemenata 
na prostorno zavisnu višegrupnu difuzionu jednačinu. Pomoću pra-
vougaonih elemenata difuziona jednačina je diskretizovana pa.je 
primenom deo po deo kontinualnlh polinoma i varijacionog metoda 
svedena na sistem linearnih algebarskih jednačina. Globalna mat-
rica je faktorizovana Choleski postupkom, a numeracija tačaka i 
konačnih elemenata linijskom disekcijom. Sumirana iskustva su is-
korišćena za pravljenje računarskog programa KONEL za rešavanje 
višegrupne dvodimenzione difuzione jednačine. 
ABSTRACT - The application of finite element method in spa-
ce dependant neutron diffusion equation is shown in this paper. 
The diffusion equation is first discretized using rectangular 
elements. Applying piecewise polynomials in a variational form 
of the diffusion equation a system of linear algebraic equati-
ons is obtained. Choleski decomposition is used to factorise 
the global matrix and numeration of nodes and finite elements is 
performed by line dissection. The multigroup two-dimensional dif-
fusion code KONEL was written based on this method. 
1. UVOD 
Primena metode konačnih elemenata na rešavanje difuzione 
jednačine se zasniva na korišćenju interpolacionih polinoma za 
prostornu aproksimaciju i varijacionog računa. Ukratko rečeno, 
to je postupak kojim se parcijalna diferencijalna jednačina dis-
kretizuje deljenjem oblasti problema u odredjen broj podoblasti, 
koje se zovu elementi. Interpolacione funkcije (obično polinomi) 
su formulisane za svaki element, po članovima koji se odnose na 
diskretne tačke. Ove tačke mogu da se nalaze u uglovima, na stra-
nama i u samoj unutrašnjosti elementa. Interpolacione funkcije 
IX .170 
susednih elemenata, odnosno sami elementi, su povezani preko nodal 
nih parametara, koji uspostavljaju uslove kontinuiteta na granica-
ma elemenata. Reaultujući deo po deo kontinualne probne funkcije 
su uključene u karakteristični funkcional problema. Postavljajući 
uslov da funkcional bude stacionaran, dobija se sistem algebarskih 
jednačina za nodalne parametre. Metod konačnih elemenata se svodi 
na Ritz-ov postupak kada su zadovoljeni potrebni uslovi kontinuite 
ta izmedju susednih elemenata. Sa istom lakoćom i tačnošću ovaj me 
tod razmatra pravilne i nepravilne geometrije, homogene i hetero-
gene materijale i bilo koju kombinaciju graničnih uslova. Zahteva-
na tačnost se dobija sa minimumom simultanih jednačina koncentrišu 
ći se na gustinu nodova, ne obraćajući pažnju na koordinatne povr-
šine u zonama gde se pojavljuju veliki gradijenti zavisnih promen-
ljivih. 
2. RAČUNARSKI PROGRAM KONEL 
U programu KONEL energetska zavisnost difuzione jednačine je 
tretirana višegrupno. Difuziona jednačina u oblasti ti sa granicom 
ЗП se piše u sledećem obliku: 
Нф = -v D уф + ru, = S. (1) 
Granični uslovi za Ф su: 
• I ,„ • 0 1 D I • 0 (2) 
1
 ЗП 
1 da je • i D — kontinualno na medjupovršinama, gde je - — 
3 n j n 
spoljašnja normala na površinu. 
Koristeći Kurganoff-ov funkcional /1/: 
1(ф) = _ j ф s d a = (Нф.ф) - 2 {S,ф) (3) 
п 
gde ( , ) predstavlja unutrašnji proizvod, dobija se 
КФ) = I (-Ф 'D + ^Ф2 - 2 S0)dC . (4) 
П 
Ako se na prvi član gornje jednačine primeni Green-ova f o r m u l a , 
jednačina dobija oblik: 
К Ф ) = (D( 'Ф^2. + 7-Ф2 - ?s$)dn- Сф i | d( jn) . (5) 3 n 
зп 
IX .171 
Prilikom inegracije zadnjeg člana gornje jednačine, a s obzirom na 
orijentaciju normale na površinu, on postaje nula /2/. 
Iz ovoga se vidi da difuzionet jednačina mora da zadovolji 
suštinske granične uslove, dok su prirodni sadržani u njoj. Ovakva 
formulacija difuzione jednačine se zove slaba formulacija. Na osno-
vu prethodnog se zaključuje da su kriterijumi za izbor bazisne funk-
cije kojom se aproksimira fluks manje restriktivni i automatski je 
izbor veći. 
Kao što je ranije rečeno ceo prostor u kome se vrši proračun 
je podeljen na konačne elemente koji su medjusobno povezani granič-
nim uslovima na njihovim dodirnim površinama. Usled toga funkcio-
nal (5) ima sada sledeći oblik: 
I(d>) = I 
I 
r, . ?
 A . 
(D (^Ф)1" + I <jr - zSfjdsi^ 
"I 
gde se I odnosi na broj konačnih elemenata. 
Pošto program KONEL vrši proračun u X-Y geometriji probna 
funkcija će imati oblik: 
Ф
 = I [ Фij wik (* ) wj t ( y )-
п 
i J 
m 
W
"
(X) =
 Да ft 1 * Hi ^ 
k;U t^j 
(7) 
Kao što se vidi za bazisnu funkciju je uzet Lagrange-ov interpola-
cioni polinom, dok je Ф^ ј vrednost fluksa u tački ij koju treba 
odrediti. 
Kada se na osnovu varijacionog principa primeni stacionar-
nost funkcionala 1(Ф) u svakoj tački ij dobija se: 
51(Ф) 
Зф. . 
IJ 
= 0 ( 8 ) 
' i = l ,n 
j = l ,m 
što predstavlja sistem linearnih algebarskih jednačina. 
Kao posledica lokalnog karaktera konačnih elemenata dobija 
se retka matrica, odnosno matrica sa puno nultih elemenata. Osim 
toga ova matrica je simetrična i pozitivno definisana. Sve ove 
osobine pružaju mogućnost da se za rešavanje sistema linearnih 
IX .172 
algebarskih jednačina primeni neka od direktnih metoda koje su znat 
no brže od iterativnog postupka. U programu je upotrebljena fakto-
rizacija Choleskl /3/. Prlmenom ove faktorizacije originalna matri-
ca prelazi u proizvod gornjotrougaone i donjotrougaone matrice. U 
obe trougaone matrice se pojavljuju novi nenulti članovi u odnosu 
na konfiguraciju originalne matrice. Broj ovih novih nenultih čla-
nova u trougaonim matricama može biti veći ili manji zavisno od to-
ga kakav način obeležavanja tačaka i konačnih elemenata se izabere. 
U programu KONEL primenjen je metod linijske disekcije prostora na 
konačne elemente (sl.l). Ovim postupkom se ne dobija najmanji broj 
novih nenultih elemenata u matricama, kao što je slučaj sa višestru 
kom disekcijom (si.2), ali sa stanovišta programiranja je znatno 
jednostavniji od drugog metoda /4/, /5/. 
Klasična Choleski dekompozicija realne simetrične pozitivno 
T definisane matrice A reda n, je data sa A = T T gde je: 
t U t 1 2 ' ' "'in t u 0 ... 
. 0 
т = o t22.. i 
T T = t 1 2t 2 2.. . 0 
0 0 . . . 
•
 -tnn 
.. t 
( 9 ) 
Množeći matrice T i T medjusobno dobijamo sledeće jednačine koje 
T: 
(i< j) 
dozvoljavaju da se odrede elementi t^j matrice 
+
 »Si^j + ••• + *И*1ј - aij 
2 2 2 
+
 <Ч2+---+ tii "ii 
(10) 
Na osnovu ovoga može se pisati sukcesivno: 
"11 11 (j >1) 
i-1 
4 i = 1/ - U - £ 
aij 
t.. * ij 
i-l 
I 
k=l 
ki 
'kl 'kj 
' u 
(1< ii n) 
(1 < ј) 
( U ) 
i >J-
IX .173 
Ako se zna da je Аф = Q, onda se mogu pisati sledeće dve jednačine: 
T T у = 0 i Тф = у (12) 
ili razvijeno u sisteme jednačina: 
4 i у1 = Q1 
t 1 2 Ух + t 2 2 У 2 = 2 ( 1 3 ) 
fcln Vi + *2п y 2 + * *пп - 9„ 
što odgovara prvoj jednačini, i 
cii ф1 • Ч г ф2 + + * т фп = 
*гг ф2 + + l2n фп - ( 1 4> 
t ф = у 
nn n J n 
što odgovara drugoj jednačini. 
Iz ova dva sistema jednačina sukcesivno se dobija: 
o, - Ti t k l У к -
1
 у. - f 'ik Фк 
ф
п
- Џ - - , \ - 1 k : 1 + 1 (i««)- (is) 
n
 nn 1 ii 
Pošto je, kao što je već rečeno matrica A simetrična a mat-
rica T gornja trougaona, potrebno je pisati samo j (n+1) gornjih 
koeficijenata а^j, i t ^ (i ^ j) kada se računa. Dakle pri prime-
ni ovog metoda, jednačine (11) i (15) se koriste za skukceslvno 
računannje koeficijenata, t ^ i у^ (i=l,2 n), a onda postup-
kom u nazad, jednačinom (16) dobi jaju se nepoznate ф^ ^ (i = n, 
n-1 1). 
Dakle, kao što se vidi kod metode konačnih elemenata vred-
nosti fluksa u tačkama se dobijaju direktnim rešavanjem sistema 
linearnih algebarskih jednačina za svaku energetsku grupu, a ne 
kao kod konačnih razlika unutrašnjim iteracijama, odnosno iteri-
ranjem po fluksu. Ovo dakako povećava brzinu rada kod metode ko-
načnih elemenata. 
IX .174 
Nalažanje sopstvene vrednosti problema, kao i kod metoda ko-
načnih razlika tako i kod metode konačnih elemenata se račune spo-
ljasnjim iteracijama, odnosno iteriranjem po izvoru. U programu 
KONEL ne koristi se nikakvo ubrzavanje spoljašnjih iteracija. 
Na početku programa direktno se učitavaju makroskopski efi-
kasni preseci i proračun se vrši u više zona. Na izlazu program da-
je efektivni faktor umnožavanja neutrona, prostornu raspodelu fluk-
sa i snage. Program je još u fazi testiranja. 
3. PRIMER PRORAČUNA I ZAKLJUČAK 
Za testiranje programa KONEL izvršen je proračun dvodimenzio-
nog IAEA BENCHMARK problema /6/, na računarskoj mašini IBM 370/3031 
Kao referentna vrednost za k ^ je uzeta vrednost dobijena ekstra— 
polacijom vrednosti koje su dobijene programom VENTURE /7/ i koja 
iznosi * 1,02959. Poredjenja radi izvršen je proračun dvodimen 
zionog IAEA BENCHMARK problema i programom VAMPIR koji se zasniva 
na primeni metode konačnih razlika. U tabelama 1. 1 2 . dati su ka-
rakteristični parametri oba programa za ovaj proračun u zavisnosti 
od reda aproksimacije odnosno koraka mreže.. Na slici 3. data je 
prostorna raspodela odstupanja srednje snage izračunate programima 
KONEL i VAMPIR od vrednosti dobijenih ekstrapolacijom rezultata iz-
računatih programom VENTURE, izražena u procentima. 
77 nrilraTarvi h i-ai-o c o mnvo a^l^ l "i ] ] ? -i i fla -ip fflpt-nH r. _ 
,— u... . ... — — ~ — - — —— vi« -~ 
načnih elemenata znatno brži od metoda konačnih razlika i da 
postiže zadovoljavajuću tačnost. 
4. LITERATURA 
/1/ A.Kavenoky, J.J.Lautard, A Finite Element Depletion Diffusion 
Calculation Method With Space Dependent Cross Sections, 
Nucl.Sci and Eng., 64, p.563, (1977) 
/2/ C.M.Kang and K.F.Hansen, Finite Element Methods for Reactor 
Analysis, Nucl.Sci. and Eng., 51, p.456 (1973) 
/3/ D.M.Davierwalla, A.Pelloni, D.Altiparmakov, Finite Elements 
In Multidimensional Nuclear Diffusion, Proc.Int.Topi.Mtg. 
Advances in Mathematical Methods for the Solution of Nuclear 
Engineering Problems,Vol.I, p.335, Munchen, 1981. 
/4/ D.Altiparmakov, Some Improvements in Finite Element Method 
Applied to Neutron Diffusion Calculations, EIR Ber.No.429,1983. 
/5/ A.George, Numerical Experiments Using-Dissection Methods to 
SolveП by n Grid Problems,SIAM J.Numer.Anal,Vol.14,No.2,1977. 
/6/ C.E.Higgs, FINELM Architecture and Usage, Newsletter of the 
NEA Data Bank, 30, p. 179 (1983) 
/7/ M.N.Zizin, L.K. Shishkov, L . IJ. Yarosl avt seva, Fizika Yaderny-
kh Reaktorov: Testovye Nejt.ronno-Fizicheskie Pasrhety 
Yademykh Reaktorov, Moskva, Atornizdat, 1980. 
• 
• 
• 
• ••
 • 
• 
• 
• 
•
п
 • 
• 
• 
п 
• 
• 
••
 • • • • • • 
•
н
 • 
• 
• 
п 
• 
• p • • • • • • • • п • п • п 
• 
• 
• 
• 
• 
• p u • • 
•
п
 • 
• 
• 
I l 
l i 
1 5 
Ш 
» • S 
ш 
i 
f 1 i 
Ш ! 
с ••••г 
г п 
L ••••п г 
г 
L UUU[]r 
г г 
L •••пг r 
г 1 ••••г 
г r 
r 
ппппг l JI 1 « 1 U_J i I г i Г " I 
U ••••Г 
г 
г 
••••C • 
с 
Si.T"XI 
IX .176 
PROGRAM KONtL 
Red 
aproksimacije Broj iteraclja Vreme proračuna 
linearna 1,0305338 80 24,19 s 
kvadratna 1,0293713 70 56,8 * 
kubna 1,0293779 95 2 *in 44,65 9 
Setvrtog rede 1,0293217 95 5 ain 19,57 « 
Tabela 1. 
PROGRAM VAMPIR 
Korak 
•ге2« 
kaff Broj iteraclja Vreae proračuna 
S ca 1,0266695 205 3 ain 52,61 s 
2,5 ea 1.0257893 178 18 ain 12,17 e 
Tabela 2. 
-28,1 
-23,5 -9,1 
-5.8 -5,7 
-3.6 
-16,6 
-6,7 -8,9 
-1,8 -5 
-3.2 
-7,5 
-6,8 0,5 
-1,5 0,8 
0 7 
-29,3 
-24,8 -9,8 
-6,2 -5,8 
5.8 
5 1,6 
1.3 0,9 
0,6 
3,8 
2,9 1,3 
0,8 1 
0.7 
-3,9 
-3,9 -0,3 
-1 0,3 
0,2 
-17,6 
—16,2 -5,4 
-4,1 -2,9 
-l.fl 
f.9 2,7 
2 1,2 
0,8 
b,9 
8,2 2 
1.5 1 
0,6 
3,2 0.8 
0,8 0,8 
0,4 
—0,4 
-2 -0,5 
-0,6 0,4 
0.2 
-6,3 
-8.4 -1,г 
-2,2 0,4 
- 0.3 
-28,7 
-25,4 -9,6 
-6,4 -5,5 
10,3 
10 3,7 
2,6 1,6 
9,2 
8,8 3,1 
2,2 1,3 
8,1 
',6 2,6 
i.9 1,2 
о.д 
6 
5.4 1,9 
1.5 1,4 
Ifl 
-0,5 -0,1 
-0,2 0.6 
-5.4 
-6,6 -1.3 
-1,8 ~ 0,2 0.1 
-19,9 
-18,4 -5,2 
-4,7 -2.4 
4,1 
6.4 -1,6 
1.5 -3.5 
-3,8 
11.9 
11.3 4,5 
3 2,1 
l.S 
10,1 
9,4 3,2 
2,4 1,4 
1 
10 
9,1 3,5 
2.4 1,6 
1 
-5,3 
0,2 -4,7 
ОД -0,6 
-3.7 
2.9 
0,7 0,7 
0,3 1,1 
-4,3 
-5,9 -1,2 
-1,6 0,3 
0.2 
-17,7 
-17,7 -5,6 
-4,6 -2.3 
prograa KONEL 
n Linearna aproksimacija 
*m kvadratna aproksinaciJa 
nm kubna aprokeiaecija 
ткм aproksimacija Četvrtog steoena-
prograa VAMPIR 
• korak lareSe 5 c» 
•* korak mreže 2,5 ci 
Slika 3. - Prostorna raspodela odstupanja srednje snage od rezultata 
dobtjenlh ekstrapolaciJom iz programa VENTURE izražena u % 
za 2 D IAEA BENCHMARK 


The application of the finite element method to the multigroup two-dimensional neutron diffusion equation on the KONEL code

http://vinar.vin.bg.ac.rs/bitstream/handle/123456789/8545/38039608.pdf?sequence=1&isAllowed=y

The application of the finite element method to the multigroup two-dimensional neutron diffusion equation on the KONEL code

Abstract

Similar works

Full text

Available Versions

VinaR - Repository of the Vinča Institute of Nuclear Sciences