Validations of the radiation transport module NEUTRO: A deterministic solver for the neutron transport equation

Valida ions agains SINBAD benchma k o he adia ion anspo module, a de e minis ic sol e o he neu on anspo equa ion. Alejand o Soba1,2, Mau icio E. Cazado2, Guillaume Houzeaux1, Albe Gu ie ez-Milla1, Me i J. Man sinen1,3 and Xa ie Saez1 1 Ba celona Supe compu ing Cen e , Ba celona, Spain. 2 CNEA-CONICET, Buenos Ai es, A gen ina. 3 ICREA, Ba celona, Spain. Keywo ds: Radia ion anspo , neu on-wall in e ac ion, usion eac o s shielding, SINBAD. Abs ac We p esen signi ican imp o emen s and alida ions o e a de e minis ic neu on anspo code (NEUTRO) dedica ed o sol ing he Bol zmann T anspo Equa ion. The code is in eg a ed as a module in he Alya sys em de eloped by he Ba celona Supe compu ing Cen e which uses he Disc e e O dina es Me hod o e an angula po ion, mul i-g oup o ene gy disc e iza ion and FEM o e uns uc u ed meshes o ea special complex domains. Ma e ial aniso opy o he sca e ing medium is in oduced in o he sca e ing ke nel using eal base exp essions o sphe ical ha monics. In o de o build he o al c oss-sec ion and he espec i e g oup ma ix o he elas ic c oss-sec ion, we use he NJOY code. We es he sol e using di e en geome ies, o de s o in eg a ion o he angula disc e iza ion and numbe o ene gy g oups. Finally, we compa e some esul s agains benchma ks ob ained om an IAEA da abase, basically measu ing leakage e ec s o se e al ma e ials. 1. In oduc ion In his pape , we p esen some imp o emen s and alida ions o a de e minis ic adia ion anspo code based on he Bol zmann T anspo Equa ion and in eg a ed in o he Alya sys em a Ba celona Supe compu ing Cen e [1], which is in ended o mul i-physics applica ions in he nuclea usion ield [2, 3]. Imp o ing he unde s anding o his ma e is essen ial o he de elopmen o u u e usion eac o s, which equi es highly demanding simula ions o he neu on anspo p ocesses. In pa icula , sol ing his kind o adia ion anspo equa ions in o de o p edic he beha iou o neu on damage and adia ion-ma e in e ac ion demands high memo y equi emen s and a high compu a ional complexi y [4]. The main idea behind his de elopmen has a double impac by in eg a ing wo complex esea ch ields. On he one hand, i will be aligned wi h he Eu opean Fusion Roadmap [5] o he ealiza ion o usion ene gy, being ocused on he i ium sel -su iciency om he Mission 4 and he design o he Tes Blanke Module (TBM). I will p o ide a cu ing-edge compu a ional ool o add essing his complex mul i-physics p oblem bo h in exis ing and u u e usion de ices which is in inc easing demand wi hin he usion communi y. On he o he hand, i has he po en ial o become a d i e o u u e usion eac o design challenges gi en i s High-Pe o mance Compu ing (HPC) capabili ies which make he analysis o di e en design a ia ions easible. This new modelling ool o usion is based on he pa allel compu a ional mechanics code Alya [6, 7]. Alya is a amewo k c ea ed o sol e di e en physical phenomena in a coupled way in la ge scale supe compu e s, accomplishing high pa allelism s anda ds and scalabili y. The code 0DQXVFULSW1HZ6XEPLVVLRQ Final e sion a ailable ia: h ps://doi.o g/10.1016/j. usengdes.2021.112497 is o ganized in a modula way including se ices, ke nel and modules. The ke nel is he co e o Alya and implemen s inpu /ou pu ope a ions, ou ines o he execu ion o he main loop and all ope a ions ela ed o he mesh geome y. Se ices cons i u e he unc ionali ies o e ed o he execu ion, o which he mos impo an is he se ice ha deals wi h he di ision o he wo k and is esponsible o he communica ion be ween nodes. Alya’s modules, as in he case o he NEUTRO module c ea ed o neu on anspo , implemen he se o Pa ial Di e en ial Equa ions (PDEs) which de ine he physics o he p oblem. The modules a e coded and compiled independen ly and can be easily coupled oge he . Using FEM in ol ing millions o billions o elemen s equi es gene al spa se sys ems o equa ions sol ed using i e a i e me hods. Mode n supe compu ing a chi ec u es consis o independen compu ing nodes, ypically mul ico e and he e ogeneous a chi ec u es, in e connec ed by a sha ed communica ion ne wo k. Alya uses he Message Passing In e ace (MPI) lib a y o implemen communica ions be ween nodes and Open Mul i-P ocessing (OpenMP) o sha ed-memo y pa allelism. De ailed desc ip ions o he s uc u e o he Alya code can be ound in e e ences [1, 6-9]. The isual pos -p ocessing o he esul s was pe o med wi h Pa a iew so wa e [10]. The sol e uses he Disc e e O dina es Me hod [11-13] o e an angula po ion, a mul i-g oup me hod o he ene gy disc e iza ion and he Fini e Elemen Me hod (FEM) [9, 14] o e uns uc u ed meshes o ea special complex domains. We use he NJOY code [15] o build he o al c oss-sec ion and he espec i e g oup ma ix o he elas ic c oss-sec ion. NJOY eads he da a om an ENDF/B da abase [16] p o ided by Los Alamos Na ional Secu i y o ice. We es he sol e in wo-dimensional geome ies (slabs o di e en dep hs), and h ee-dimensional geome ies using laye s and sphe es. Quad a u es o se e al o de s in he equal-le el disc e e o dina e se a e a emp ed [17, 18] and due o he sca e ing o high ene gy 14 MeV usion neu ons being highly aniso opic, a new ype o quad a u e based in he FEM disc e iza ion [19] is conside ed. In addi ion, we ha e es ed he sol e using di e en ene gy g oups disc e iza ions, om 3 o 256, and he esul s ob ained we e alida ed agains o he codes and he SINBAD (Shielding In eg al Benchma k A chi e Da abase) [20] benchma ks, ob ained om an IAEA da abase, in which he leakage e ec s o se e al ma e ials we e measu ed. Al hough he alida ions included in his wo k a e compa ed wi h he neu on anspo esul s, he code is p epa ed o sol e p oblems o gamma adia ion anspo as well [21, 22]. Valida ions o he la e will be p esen ed in u u e wo ks. I is impo an o men ion ha al hough he e a e many codes in he wo ld ha sol e he adia ion anspo equa ion [23-26], none o hem is in eg a ed wi h compu a ional mechanics amewo ks o mul i-physics coupled p oblems and hey ha e been de eloped as s and-alone so wa e. In his sense, he in eg a ion o he NEUTRO module in o he Alya sys em has an inno a i e aspec wi h espec o o he be e -known solu ions. The inal in en ion o his de elopmen is o achie e a mul i-physics ea men o high- olume en i onmen s composed o mul iple ma e ials, being able o sol e hei he mal, mechanical and neu onic beha iou a he same ime, using he coupling capabili ies p o ided by he Alya sys em. In his espec , we made a deep scalabili y s udy o he sol e , showing he beha iou o he code and he quali y o he esul s as we inc ease he numbe o p ocesso s in ol ed in he p oblem sol ing. Some o hese scalabili y measu es a e included in sec ion 3. This pape is o ganized as ollows: in sec ion 2, we b ie ly desc ibe he sol ing amewo k o he NEUTRO module included in Alya. In sec ion 3, we include some o he es s ca ied ou wi h his module and we show some nume ical d awbacks ha had o be sol ed in o de o achie e a good quali y in he inal esul s. In addi ion, we p esen some scalabili y esul s including he beha iou o he code when a la ge numbe o p ocesso s a e in ol ed in he calcula ions. In sec ion 4, we show he compa isons agains some selec ed expe imen s o he SINBAD da abase. Finally, in sec ion 5, we p esen ou conclusions. 2. Nume ical adia ion anspo equa ion The gene al s a iona y adia ion anspo equa ion sol ed a he momen in he NEUTRO module included in Alya has he ollowing exp ession [7-10]: Ω∙∇φ+Σ$(%,&,Ω)'=* * Σ-(%,&.→&,Ω.→Ω)'0Ω.0&. 12 3 4 (1) Whe e Ω: Solid angle o he possible di ec ion in whe e a pa icle a els. E: Pa icle ene gy. : Gene al spa ial a iable. φ: Flux as a unc ion o , E and Ω. Σ$( ,E,Ω): To al mac oscopic c oss-sec ion. Σ-( ,E.→E,Ω): To al mac oscopic sca e ing c oss-sec ion. The ene gy mul i-g oup app oxima ion [9, 27, 28] and he disc e iza ion o angula dependence we e applied o equa ion (1). Fo he la e , he Me hod o Disc e e O dina es (equal-le el se , Sn) [11-14] o quad a u es se based on discon inuous ini e elemen s [19] may be chosen. The mul i-g oup me hodology consis s o di iding he neu on ene gy ange in se e al g oups (G) and e alua e all he ene gy-dependen pa ame e s in ela ion o hose ene gy g oups G. Disc e iza ion o e he angula (di ec ional) a iable is usually pe o med by means o mechanical quad a u e. The con inuous angula a iable Ω is eplaced by a se o disc e e di ec ions, each o which has a ela i e weigh . The di ec ions may be ep esen ed as poin s on he uni sphe e. Using a mul i-g oup app oxima ion and an angula disc e iza ion, he in eg al-di e en ial anspo equa ion is con e ed in o a sys em o G×M di e en ial equa ion as ollows: 78∙∇φ98 +Σ$;%,&9<φ98 = > > Σ-9?,9 8?,8(%)φ9?8? @ 8?AB 8.C8 D 9?AB 9.C9 (2) Wi h F =1,…,I; K=1,…,L Ma e ial aniso opy o he sca e ing medium is in oduced in o he sca e ing ke nel using eal- based exp essions o sphe ical ha monics. Combining his app oxima ion wi h he ene gy g oup di ision and disc e e o dina e me hods, he ollowing exp ession o he sca e ing e m is ob ained:     4 2!+1Y"#"(Ω#) " #"%&" ' "%* , -.%/ 0 3.%/ Y 5"#"(Ω#6) Σ#"3.,3 -6,-(:)φ36-6(:) (3) ( he o de o he sphe ical ha monics was changed o ‘ml’ o a oid con usion). In his exp ession, a maximum alue o L is chosen o cu he expansion. I L=0, he ke nel exp ession in equa ion (2) is eco e ed. Any alue o L g ea e han ze o in oduces angula aniso opy. The mos gene al bounda y condi ion o he anspo equa ion is he Albedo condi ion [1,9], which consis s o : φ<=( ,E,Ω!")= φ!" + αφ( , E, Ω%) (4) whe e φ!" is a ixed lux, α akes he alue 0 o he acuum bounda y condi ion and α=1 o he gene al e lexi e condi ion. Gene ally, o model a e lec o su ounding he co e he alue o α is 1. O e all, he Albedo coe icien α may a y wi h he e lec o p ope ies and should also cap u e angula and ene ge ic edis ibu ion o he e lec ed neu ons due o hei di usion h ough he e lec o . The e a e wo main s eps when implemen ing mi o bounda y condi ions. Fi s ly, i is necessa y o calcula e he no mal ec o o he su ace on o which he ays a e e lec ed. Secondly, he calcula ion o he e lec ion angle is equi ed. Since we a e wo king wi h disc e e o dina es we need o p ojec he ajec o y on o one o he disc e ized di ec ions ha we ha e. In each bounda y whe e he e lec ing condi ion is u ned on, he algo i hm de eloped de e mines which di ec ion is he closes o he cu en di ec ion o e lec ion o use i in he building o he bounda y con ibu ion. In he disc e e o dina es me hod, gi en a di ec ion used in he angula disc e iza ion, i is necessa y o e alua e he co esponding specula e lec i e di ec ion a e he neu on a i ed a he e lec i e medium. To do so, a special unc ion was p og ammed in o de o de e mine a he beginning o he simula ion he co esponding di ec ion o each no mal p esen in he limi s o he domain o ou p oblem. Using his p ocedu e, bo h angula di ec ion and e lec i e di ec ion a e e alua ed o e each bounda y whe e a e lec i e condi ion is equi ed. In o de o build he o al, cap u e and elas ic sca e ing c oss-sec ions and he espec i e g oup ma ix o elas ic and ission c oss-sec ion, he NJOY code [15] was used. NJOY eads and p ocesses da a om ENDF/B da abase [16]. F om SINBAD benchma k da abase, some exe cises des ined o analyse he leakage e ec o se e al ma e ials we e selec ed. Ini ially, he inpu iles wi h he o al and sca e ing ma ix o i on, nickel, aluminium, ungs en, manganese and silicon wi h NE-ene gy g oups we e c ea ed. We es ed wi h NE = 238, 118 and 59 g oups o compa e wi h he benchma k cases. Fo mos o he shielding expe imen s aken om he SINBAD da abase, we used 59 ene gy g oups. In addi ion, o es ing coa se de ini ions, 29- ene gy g oups and 3-ene gy g oups we e used. 3. Nume ical d awbacks and pe o mance analysis 3.1 Angula disc e iza ion: Once he de elopmen o he mul i-g oup equa ion was comple ed and be o e going o he alida ion s age, we ca ied ou a se ies o nume ical expe imen s o e i y he unc ioning o his new e sion o he code wi h espec o wha was p e iously p og ammed, as well as i s in e ac ion wi h he Alya sys em. In he i s place, he solu ion o a sphe e wi h an in e nal hole wi h an inne lux in he adial di ec ion and acuum in he ex e nal bounda y was deeply analysed. We used a e ahed al mesh o linea elemen s o disc e ize he eal space. As he in en ion was o analyse he solu ion wi h he di e en quad a u e o mulas implemen ed in he code, a h ee-ene gy g oups disc e iza ion was used. Fo compa ison, we analysed he solu ion wi h di e en app oxima ions o he equal le el me hods and he implemen ed angula disc e iza ion based on he ini e elemen me hod wi h 512 poin s (called in e nally in he code FE64). This me hod consis ed o he disc e iza ion o he medium plane o he sphe e using a iangula mesh and he p ojec ion o e he sphe e su ace o each iangle wi h a well-de ined weigh [19]. In Figu e 1, we p esen a gene al iew o he code capaci ies. Figu e 1 a) shows a domain di ision o be dis ibu ed among 16 nodes, used o sol e he eal space using FEM. Figu e 1 b) p esen s he neu onic cu en in he z di ec ion showing he symme y o he solu ion ob ained. Figu e 1 c) is a con ou zone o one ene gy g oup as an example solu ion o he domain. All hese solu ions we e ob ained using he mo e accu a e angula disc e iza ion (FE64). Las ly, in Figu e 1 d), we plo he o al lux in all he nodes on he su ace o he sphe e o he di e en angula disc e iza ions implemen ed in he code. I can be app ecia ed how he solu ions ob ained wi h he equal le el me hod Sn4 (24 poin s) and Sn8 (80 poin s) p esen ed a g ea dispe sion and nega i e alues o he solu ion. When we passed o Sn10 (120 poin s) and Sn12 (168 poin s) posi i e alues and less dispe sion in he solu ion we e ob ained. Finally, o he quad a u e based on FE he esul s achie ed show he leas dispe sion, ob aining he uni o mi y ha he symme y o he solu ion p oposes. Le us highligh ha he Sn12 me hod has 168 in eg a ion poin s in all he sphe e while he FE64 has 512 in eg a ion poin s, inc easing he compu a ion ime o he p oblem p opo ionally. Bu he mo e impo an esul is ha using he highes quad a u e solu ion, we a oid he spu ious nega i e alues in he lux solu ion, a ecu ing p oblem in he use o disc e e me hods o sol e he anspo equa ion. Figu e 1: a) Pa i ioned mesh o pa allelize he solu ion o he eal space. b) Neu onic cu en in z di ec ion. c) Flux solu ion o one selec ed ene gy g oup d) Dispe sion o he su ace solu ion o di e en modes o angula disc e iza ion. 3.2 Compa ison be ween 2D and 3D solu ions: An impo an and e ec i e way o accele a e a p oblem om he nume ical poin o iew is o educe he dimensions o he domain o be analysed when some ypes o symme y a e p esen . This educ ion, in he case o anspo p oblems, in ol es educing he physical domain bu also he angula disc e iza ion by omi ing ou oc an s o angula di ec ions. On he one hand, his sys em inc eases he esolu ion speed, bu a he same ime enla ges he e o o he solu ion. In his sub-sec ion, we compa e a solu ion ob ained on a h ee-dimensional laye and a slab o equal wid h, bu wo-dimensional. Bo h cons i u ed by he same ma e ial, and wi h he same inne lux om he adial di ec ion (x in he slab2d and z in he 3D laye ). A acuum bounda y condi ion is conside ed in he ou pu di ec ion and e lec i e bounda ies in he es o he limi s. Figu e 3 shows he o al lux ob ained on he ex e nal ace as a unc ion o one o he signi ican coo dina es. The solu ion is di ided by he maximum ob ained in he h ee-dimensional case, eaching a di e ence o 14% in he esul be ween bo h simula ions. Figu e 2: Solu ion o a wo-dimensional slab compa ed wi h a h ee-dimensional laye o he same wid h 3.3 Pe o mance Analysis: Be o e discussing he pe o mance analysis, a ew ema ks should be made abou he pa alleliza ion in Alya. Fi s o all, he implemen ed pa alleliza ion s a egy is he main-seconda y echnique. This means ha when he pa alleliza ion akes place, he main p ocess is in cha ge o dis ibu ing he da a and coo dina ing he o he p ocesses (seconda y). The seconda ies a e he ones ca ying ou he calcula ions. All p ocesses communica e among hemsel es h ough messages using he MPI s anda d, which p o ides synch oniza ion and communica ion unc ionali ies be ween a se o p ocesses wi h dis ibu ed memo y (meaning ha no all p ocesses can access he same memo y space and he e o e need o exchange da a among hem). The Alya sys em is w i en in FORTRAN language: MPI is implemen ed as a lib a y wi hin he code om which we can call he unc ions we may need once we ha e included i , as we would do wi h any o he lib a y. Ano he impo an ea u e o he Alya implemen a ion is he possibili y o mul iplying he ini ial mesh. Using his op ion each ini e elemen is di ided in o up o 8 ini e elemen s o he same cha ac e is ics: his implies an inc ease in he numbe o elemen s o he inpu mesh. In o de o es he pe o mance, a h ee-dimensional sphe e o Fe56 ma e ial wi h an in e io sou ce in he adial di ec ion using 59-ene gy g oups and Sn12 o he angula quad a u e (see igu e 1) we e used. Mul iplying by 2 and by 3, meshes o 77984 and 623872 ini e elemen s a e de ined, espec i ely. In Figu e 3 he speed up o each case is p esen ed. In he case o 77984 elemen s, a 95% e iciency wi h 256 p ocesso s is achie ed, while in he case o 623872 elemen s, an 82% e iciency wi h 1024 p ocesso s and 70% wi h 2048 p ocesso s a e ob ained. Bo h cases show sa u a ion due o he inc ease in communica ion. No e ha his e iciency is measu ed wi h espec o he un on 16 ( o he small case) and 128 ( o he big case) MPI p ocesses. Despi e he sa u a ion men ioned, he esul s show good pe o mance wi h mo e han wo housand p ocesso s. Fo all he calculus we used a gene alized minimal esidual algo i hm (GMRES) [29] o sol e he nonsymme ic linea sys em p oduced by his ype o p oblem. We used a diagonal p econdi ione , which esul ed in he mos sui able op ion o each con e gence in less ime. Figu e 3: Speedup o wo cases, using wo meshes, 77984 elemen s wi h an a e age olume o 6.73x10-6 and 623872 elemen s wi h 8.4x10-7 a e age olume. 4. SINBAD es ing In o de o alida e he de eloped NEUTRO module, we selec ed a numbe o benchma k cases o es ing pu poses om SINBAD, an OECD/NEA Da a Bank (O ganiza ion o Economic Coope a ion and De elopmen , Nuclea Ene gy Agency Da a Bank) and ORNL/RSICC (Oak Ridge Na ional Labo a o y, Radia ion Sa e y In o ma ion Compu a ional Cen e) join p ojec [20, 30]. The main objec i e o hese es s was o p oduce a da abase con aining an in e na ionally es ablished se o adia ion shielding and dosime y da a ela i e o expe imen s ele an o eac o shielding, usion blanke neu onics and accele a o shielding, o e ing a s anda dized ool o alida ion o adia ion anspo compu e codes and nuclea da a lib a ies. Ou pa icula analysis was based on he se o expe imen s ha ook place in Osaka Uni e si y a he OKTAVIAN acili y. De ailed desc ip ions o he expe imen s, he de ec o sys em and he sou ce spec um used can be ound in e e ences [31-36]. In summa y, his g oup o es s used hollow sphe es wi h an inne neu on sou ce wi h a p onounced peak in he ange o 14 MeV. The sphe es we e made o di e en ma e ials (Fe, Ni, Al, W, Mn, Si, among o he s) and wi h di e en diame e s in a ange o 14 cm o 100 cm. The main pu pose o hese benchma ks was o s udy he neu on and gamma- ay leakage spec a om he ou e su ace o he sphe e. The e a e wo kinds o measu emen s ca ied ou depending on he posi ion o he de ec o . Fo some ma e ials (Mn, Si, Al, W), he neu on de ec o was loca ed wi h a 55° angle wi h espec o he cen al axis o he sphe e [35, 36]. In o he cases, such as o Fe and Ni, he de ec o did no co e he en i e su ace o he sphe e bu he solid angle o 17.28° om he sphe e cen e due o he p esence o he collima o s [32-34]. In o de o pe o m he es , wo ypes o models we e used: a wo-dimensional model, using slabs wi h a wid h equal o he sphe ical shell hickness; a h ee-dimensional laye wi h a wid h equal o he sphe ical shell hickness and a h ee-dimensional sphe e wi h an in e nal hole. Table I p esen s a summa y o he es s pe o med and he ma e ials used. Table I: Desc ip ion o he Benchma k used in his wo k [31,36] Model Ma e ial Ex . Radi us (cm) In e nal hole adius (cm) E ec i e hickness (cm) De ec o posi ion Sphe e o 40 cm Silicon 20 10 10 55º de ec o Sphe e o 60 cm Silicon 30 10 20 55º de ec o Sphe e o 60 cm Mangane se 30 10 20 55º de ec o Sphe e o 40 cm Aluminium 20 10 10 55º de ec o Sphe e o 40 cm Tungs en 20 10 10 55º de ec o Sphe e o 32 cm Nickel 16 2 14 17.28º solid angle om he sphe e cen e Sphe e o 100 cm I on 50 10 40 17.28º solid angle in he sphe e cen e In Figu e 4, compa isons be ween expe imen al esul s and se e al benchma ks un wi h Alya using wo-dimensional slabs a e p esen ed. In each case, he wid h o he slabs ep esen s he e ec i e hickness o he sphe e in each expe imen . Al hough o he cases o aluminium o 10 cm hickness (Figu e 4 a), manganese o 20 cm hickness (Figu e 4 b), i on o 40 cm hickness (Figu e 4 c) and ungs en o 10 cm hickness (Figu e 4 d) he esul s show some dispe sion in he ange o ene gy be ween 0.01-1 MeV, he solu ion ob ained using he NEUTRO module shows a good ep esen a ion o he expe imen al da a endency. Figu e 4: Expe imen al measu ed and Alya esul s o Leakage s. Ene gy o SINBAD es using wo-dimensional slabs. Be ween pa en heses he wid h o he slabs conside ed, co esponding o he hickness o he sphe e on which he expe imen s we e ca ied ou . 59 ene gy g oups we e used o ene gy disc e iza ion. Figu e 5 p esen s he compa ison be ween expe imen al da a and calcula ed alues o h ee- dimensional domains, bo h laye s and sphe es. In each case, we used meshes o linea hexahed a o he laye s and linea e ahed a o he sphe es. The Sn8 quad a u e o he laye and FE64 o he sphe es we e used. 59 ene gy g oups we e used o he ene gy disc e iza ion. The esul s ob ained o a silicon laye o 20 cm hickness (Figu e 5 a), a manganese laye o 20 cm hickness (Figu e 5 b), a ungs en sphe ical shell o 20 cm hickness (Figu e 5 c) and a nickel sphe ical shell o 14 cm hickness (Figu e 5 d) show some dispe sion in di e en anges o ene gy, depending on he expe imen s. In gene al, he esul s o h ee-dimensional geome ies a e hose ha p esen g ea e dispe sion wi h espec o he expe imen al measu emen s. And among hese, hose co esponding o sphe ical geome ies a e he a hes om he measu ed alues. The e a e se e al easons ha con e ge o hese esul s. The i s is he e o ha an uns uc u ed mesh on a ci cula geome y in oduces, which in luences he en i e calcula ion om he in e nal bounda y condi ion, whe e he incoming low is ixed, o he ou e edge o he domain. Secondly, he ex ension o he domain in ol ed, which, al hough i is ea ed by a la ge numbe o elemen s, is s ill mo e di icul o ep esen by means o a mesh o e eal space. This di icul y is pa ly a esul o using wo-dimensional domains, which allow dense meshing wi hou aising he compu a ional cos . Fu he mo e, as in he laye s, he slabs do no ha e a meshing e o and i is possible o co ec ly ep esen he en i e analysed su ace.