APRIL 2021 www.fe.um.si/en/jet.html 2 JET VOLUME 14 / Issue 1 Revija Journal of Energy Technology (JET) je indeksirana v bazah INSPEC© in Proquest’s Technology Research Database. The Journal of Energy Technology (JET) is indexed and abstracted in database INSPEC© and Pro­quest’s Technology Research Database. 4 JET JOURNAL OF ENERGY TECHNOLOGY Ustanovitelj / FOUNDER Fakulteta za energetiko, UNIVERZA V MARIBORU / FACULTY OF ENERGY TECHNOLOGY, UNIVERSITY OF MARIBOR Izdajatelj / PUBLISHER Fakulteta za energetiko, UNIVERZA V MARIBORU / FACULTY OF ENERGY TECHNOLOGY, UNIVERSITY OF MARIBOR Glavni in odgovorni urednik / EDITOR-IN-CHIEF Jurij AVSEC Souredniki / CO-EDITORS Bruno CVIKL Miralem HADŽISELIMOVIC Gorazd HREN Zdravko PRAUNSEIS Sebastijan SEME Bojan ŠTUMBERGER Janez USENIK Peter VIRTIC Ivan ŽAGAR Uredniški odbor / EDITORIAL BOARD Dr. Anton BERGANT, Litostroj Power d.d., Slovenia Izr. prof. dr. Marinko BARUKCIC, Josip Juraj Strossmayer University of Osijek, Croatia Prof. dr. Goga CVETKOVSKI, Ss. Cyril and Methodius University in Skopje, Macedonia Prof. dr. Nenad CVETKOVIC, University of Nis, Serbia Prof. ddr. Denis ĐONLAGIC, University of Maribor, Slovenia Doc. dr. Brigita FERCEC, University of Maribor, Slovenia Prof. dr. Željko HEDERIC, Josip Juraj Strossmayer University of Osijek, Croatia Prof. dr. Marko JESENIK, University of Maribor, Slovenia Izr. prof. dr. Ivan Aleksander KODELI, Jožef Stefan Institute, Slovenia Izr. prof. dr. Rebeka KOVACIC LUKMAN, University of Maribor, Slovenia Prof. dr. Milan MARCIC, University of Maribor, Slovenia Prof. dr. Igor MEDVED, Slovak University of Technology in Bratislava, Slovakia Izr. prof. dr. Matej MENCINGER, University of Maribor, Slovenia Prof. dr. Greg NATERER, Memorial University of Newfoundland, Canada Prof. dr. Enrico NOBILE, University of Trieste, Italia Prof. dr. Urška LAVRENCIC ŠTANGAR, University of Ljubljana, Slovenia Izr. prof. dr. Luka SNOJ, Jožef Stefan Institute, Slovenia Izr. prof. dr. Simon ŠPACAPAN, University of Maribor, Slovenia Prof. dr. Gorazd ŠTUMBERGER, University of Maribor, Slovenia Prof. dr. Anton TRNIK, Constantine the Philosopher University in Nitra, Slovakia Prof. dr. Zdravko VIRAG, University of Zagreb, Croatia Prof. dr. Mykhailo ZAGIRNYAK, Kremenchuk Mykhailo Ostrohradskyi National University, Ukraine Prof. dr. Marija ŽIVIC, University of Slavonski Brod, Croatia 6 JET Tehnicni urednik / TECHNICAL EDITOR Sonja Novak Tehnicna podpora / TECHNICAL SUPPORT Tamara BRECKO BOGOVCIC Izhajanje revije / PUBLISHING Revija izhaja štirikrat letno v nakladi 100 izvodov. Clanki so dostopni na spletni strani revije -www.fe.um.si/si/jet.html / The journal is published four times a year. Articles are available at the journal’s home page - www.fe.um.si/en/jet.html. Cena posameznega izvoda revije (brez DDV) / Price per issue (VAT not included in price): 50,00 EUR Informacije o narocninah / Subscription information: http://www.fe.um.si/en/jet/ subscriptions.html Lektoriranje / LANGUAGE EDITING Terry T. JACKSON Oblikovanje in tisk / DESIGN AND PRINT Fotografika, Boštjan Colaric s.p. Naslovna fotografija / COVER PHOTOGRAPH Jurij AVSEC Oblikovanje znaka revije / JOURNAL AND LOGO DESIGN Andrej PREDIN Ustanovni urednik / FOUNDING EDITOR Andrej PREDIN Izdajanje revije JET financno podpira Javna agencija za raziskovalno dejavnost Republike Slovenije iz sredstev državnega proracuna iz naslova razpisa za sofinanciranje domacih znanstvenih periodicnih publikacij / The Journal of Energy Technology is co-financed by the Slovenian Research Agency. Spoštovani bralci revije Journal of energy technology (JET) Energija vetra je posledica obsevanja Zemljine površine s soncnimi žarki. Potenciali vetrne energije v svetu so zelo veliki, mnogo vecji kot trenutno izkorišcanje energije vetra. Intenzivnost vetrne energije je geografsko razlicna, pravzaprav pa ni države oz. pokrajine, kjer ni obmocij s potenciali za izkorišcanje. Tovrstno energijo so znali izkorišcati že v preteklosti; tako naj bi že 5000 let pred našim štetjem uporabljali vetrno energijo za pogon ladij vzdolž reke Nil, pred približno 200 leti pred našim štetjem pa so jo Kitajci uporabljali za pogon enostavnih crpalk. Uporaba energije vetra se je nato razširila. Tako so ljudje na ozemlju srednjega vzhoda v 11. stoletju že uporabljali energijo vetra za pogon vetrnih crpalk in vetrnih mlinov. Še posebej velikopotezna je bila uporaba vetrne energije na Nizozemskem za izsuševanje potiskanja morja. Danes predstavlja vetrna energija v svetu zelo pomemben clen v proizvodnji elektricne energije. V letu 2020 imajo inštalirane vetrne elektrarne moc že preko 750 GW. Globalno najvec energije iz­korišcajo na Kitajskem, v ZDA in Nemciji. Najvec elektricne energije na prebivalca pridobijo na Dan-skem, Švedskem, Irskem in v Nemciji. Najvecje vetrne elektrarne imajo moc okoli 10 MW, s premeri rotorjev okoli 160 m. V Sloveniji je izraba vetrne energije za proizvodnjo elektricne energije zelo majhna, prav zato so raziskave na tem podrocju izjemnega pomena. Vsekakor je potrebno v bližnji prihodnosti v energetiki dati vecji poudarek izrabi vetrne energije. V tej številki revije je predstavljen clanek na to tematiko, zato vas vabim k branju revije. Jurij AVSEC odgovorni urednik revije JET 8 JET Dear Readers of the Journal of Energy Technology (JET) Wind energy results from the irradiation of the Earth's surface with solar rays. The global potential of wind energy is massive, much greater than the current exploitation of it. The intensity of wind energy is geographically different, but there is no country or province in the world without areas with interesting potentials for wind energy. People have known how to harness the power of the wind since ancient times; Egyptians used wind energy to power ships along the Nile 5,000 years BC, and about 200 BC the Chinese used wind energy to power simple pumps. The use of wind energy then spread worldwide; in the 11th century, people in the Middle East used wind energy to power wind pumps and windmills. The use of wind energy in the Netherlands to recover land from the sea was particularly ambitious. Today, wind energy is an essential link in global electricity production. In 2020, the installed wind power plants had a capacity of over 750 GW. China, the United States and Germany produce the most energy. Denmark, Sweden, Ireland and Germany have the most electricity per capita from wind energy. The largest wind farms have a capacity of about 10 MW, with rotor diameters of about 160m. In Slovenia, the use of wind energy for electricity production is very small, which is why re­search in this field is extremely important. In any case, greater emphasis needs to be placed on the use of wind energy in the near future. This issue of the magazine presents an article on this very topic, so I invite readers to read further. Jurij AVSEC Editor-in-chief of JET Table of Contents / Kazalo New technique to evaluate the overall heat loss coefficient for a flat plate solar collector Nova tehnika za oceno celotnega koeficienta toplotne izgube za plošcati soncni kolektor Amor Bouhdjar, Hakim Semai, Aissa Amari. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .11 Towards Forty Years of Krško NPP Operation – An Overview of Population Exposures to Radiation Ob koncu cetrtega desetletja delovanja ne Krško - pregled sevalne izpostavljenosti prebivalstva Matjaž Koželj . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .27 Comparison of cavitation models for the prediction of cavitation around a hydrofoil Primerjava kavitacijskih modelov za numericno napoved kavitacije na hidrodinamicnem profilu Marko Pezdevšek, Ignacijo Biluš, Gorazd Hren . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .41 Modelling of magnetic regenerator and heat transfer agent microchannels Modeliranje magnetnega regeneratorja in sredstva za prenos toplote v mikrokanalih Dorin Botoc, Bianca-Eliza Oneata, Ionut Rusu, Alexandru Salceanu, Jurij Avsec. . . . . . . . . . . . .57 Instructions for authors . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .69 10 JET JET Volume 14 (2021) p.p. 11-25 Issue 1, April 2021 Type of article 1.01 www.fe.um.si/en/jet.html NEW TECHNIQUE TO EVALUATE THE OVERALL HEAT LOSS COEFFICIENT FOR A FLAT PLATE SOLAR COLLECTOR NOVA TEHNIKA ZA OCENO CELOTNEGA KOEFICIENTA TOPLOTNE IZGUBE ZA PLOŠCATI SONCNI KOLEKTOR Amor Bouhdjar1,R, Hakim Semai1, Aissa Amari1 Keywords: flat plate solar collector, overall heat loss coefficient, heat removal factor Abstract Low-temperature solar systems mostly use flat plate solar collectors. Good design and correct di­mensioning of a solar heat generator are based on precise knowledge of the characteristics of the flat plate solar collector on site. The present work considers a flat plate solar air collector, a flat plate solar water collector, and a flat plate solar water collector with air absorber cooling. The investiga­tion intends to shed light on a procedure to determine the overall heat loss coefficient and the heat removal factor using recorded system temperatures, operating parameters, and environmental data. The Hottel-Whillier-Bliss equation gives the collector useful energy. This expression is used to gener­ate a correlation for the collector efficiency through a linear fitting. We calculate the overall heat loss coefficient of the collector from the slope of the collector efficiency curve. However, we need to know the heat removal factor of the collector. In this study, we present a new technique to calculate the heat removal factor. Then we deduce the collector overall heat loss coefficient. Results show that, very often, the overall heat loss coefficient for the flat plate solar air collector and the flat plate solar water collector determined with this new method is higher than the one calculated with the empirical formula proposed by Klein. R Corresponding author: Dr Amor Bouhdjar, Tel.: +231 (0) 23189051, mailing address: PB. 62 Route de l’Observatoire Bouzareah 16340 Algiers Algeria, E-mail address: a_bouhdjar@yahoo.com; a.bouhdjar@cder.dz 1Renewable Energy Development Centre, Solar thermal and Geothermal Energy Division, Route de l’Observatoire Bouzareah 16340 Algiers Algeria However, the experimental overall heat loss coefficient for the flat plate solar water collector with air absorber cooling is smaller than the one calculated with the empirical formula proposed by Klein. The analysis shows that the overall heat loss coefficient determined with the new technique seems more realistic since all phenomena occurring during the heat transfer from solar irradiance incident on the absorber plate and transmitted to the transport fluid are considered. Povzetek Nizkotemperaturni soncni sistemi vecinoma uporabljajo soncne kolektorje z ravno plošco. Dobra zasnova in pravilne dimenzije solarnega generatorja toplote temeljijo na natancnem poznavanju znacilnosti plošcatega soncnega kolektorja na lokaciji. Pricujoce delo obravnava ravno zracni plošcati soncni kolektor, vodni plošcati soncni kolektor in plošcasti soncni kolektor vode z zracnim absorberjem. Analiza v clanku namerava osvetliti postopek dolocanja celotnega koeficienta toplotnih izgub in faktorja odvajanja toplote z uporabo zabeleženih temperatur sistema, obratovalnih parametrov in okoljskih podatkov. Enacba Hottel-Whillier-Bliss daje raziskovalcem koristne informacije. Ta izraz se uporablja za ustvarjanje korelacije za ucinkovitost kolektorja z linearno vgradnjo. Skupni koeficient toplotnih izgub kolektorja izracunamo iz naklona krivulje izkoristka kolektorja. Poznati pa moramo faktor odvajanja toplote kolektorja. V tej študiji predstavljamo novo tehniko za izracun faktorja odvajanja toplote. Nato izracunamo koeficient celotne toplotne izgube kolektorja. Rezultati kažejo, da je skupni koeficient toplotne izgube za ravno plošcasti zracni soncni kolektor in vodni plošcati soncni kolektor, dolocen s to novo metodo, vecji od tistega, izracunanega z empiricno formulo, ki jo je predlagal Klein. Vendar je eksperimentalni skupni koeficient toplotne izgube za vodni plošcni soncni kolektor z zracnim absorberjem manjši od tistega, izracunanega z empiricno formulo, ki jo je predlagal Klein. Analiza kaže, da se zdi skupni koeficient toplotne izgube, dolocen z novo enacbo, bolj realen, saj so upoštevani vsi pojavi, ki se pojavijo med prenosom toplote zaradi soncnega obsevanja, ki vpade na plošco absorberja in se prenese v transportno tekocino. INTRODUCTION With the trend of decarbonization, many countries are moving toward using alternative energy sources such as solar, wind, biomass, geothermal, and others. For the supply of heat, solar water heating systems (SWHS) are widely used, especially in countries with good solar radiation. These systems are sufficiently mature to replace gas and other fossil fuel to supply heat in domestic or industrial applications. Solar air collectors are also integrated into many heating systems such as air space heating, drying applications, such as agricultural products, timber, biomass cultivation, waste biomass, building materials and desalination, and regulating microclimate in agricultural products storage facilities. Improvement of system efficiency and managerial flexibility of the energy obtained can expand the use of solar energy. In a scenario still characterized by strong growth in the installed solar collector capacity, [1], even relatively small improvements may lead to a large increase in the overall energy production in absolute terms. Rigorous studies may also lead to adequate previsions. For this reason, research and development into the optimization of the collector characteristics and solar field design play key roles. Flat plate solar collectors (FPSC) are the main and most used solar water heating systems component. It is made of an absorber plate covered above by a transparent sheet and surrounded by an isolating material. Incident solar radiations onto the absorber are absorbed by the latter and transferred as heat to a fluid flowing within the absorber. Therefore, a flat plate solar water collector is a special type of heat exchanger. In contrast, a flat plate solar air collector (FPSAC) is a simple solar heating system with a low convective heat transfer coefficient between the absorber surface and the flowing air, resulting in a low heat transfer rate. Low thermal conductivity of air and high heat loss to the environment are the drawbacks of solar air collectors, [2, 3]. Solar thermal flat plate collector performance depends very much on the absorbed solar radiation and the energy lost to the environment. Many authors have examined various configurations to assess the influence of geometrical characteristics and operating conditions on the thermal efficiency of flat plate solar water collectors. Important parameters such as absorber thickness, riser position, the shape of the tube’s cross-section, plate material, coating effect on absorptivity, cover transmissivity, fluid properties, and mass flow rate were investigated [4]. Other studies investigate cover absorber distance, heat transfer between absorber plate and cover sheet, the roughness of absorber plate surface, or use of transparent insulation material, [5]. Lowering heat losses from the absorber to the surrounding environment is an important issue for FPSC. The use of selective coating techniques was an innovation to reduce thermal radiation heat transfer from the absorber while maintaining high plate absorptivity, [6]. Others investigated heat convection between the glass cover and absorber plate and the influence of the distance between these two components, [7]. Several other studies have been made on flat plate solar water heaters, [8-11], mainly on the determination of the global heat loss coefficient and the heat removal factor of the collector. Klein, [12], proposed a mathematical model, based on the theory of Bliss-Whillier, [13, 14], which estimates the collector efficiency factor and the heat removal factor of the collector considering many simplifying assumptions. The performance of solar air collectors is usually affected by the low heat transfer coefficient between the air and the absorber plate, [10]. The choice of the optimal collector depends on the temperature level required by the specific application and on the climatic conditions at the installation site. Therefore, in terms of efficiency, each collector displays features, making it most suitable to a given application. Several new applications of solar energy have appeared and intensified their use, as shown by works developed by Esen, [15, 16]. Efforts have been made to combine a number of the most important factors into a single formulation to have a mathematical model, which will describe the collector’s thermal performance in a computationally efficient manner. Evaluating the thermal loss coefficients is the fundamental task to assess the flat plate solar collector performance. The FPSC global heat loss coefficient UC (W/m2.K) is the sum of the top loss (Ut), the bottom loss (Ub), and the edge loss (Ue) coefficients. This loss value is established between the collector and its surrounding by conduction, infrared solar radiation, and convection heat transfer. A good approximation of the FPSC global heat loss coefficient will lead to an effective solar water system design. The present work aims to determine the overall loss heat coefficient of the flat plate solar collector using a new technique to evaluate the heat removal factor of the collector, taking into consideration experimental data to minimize the effect of assumptions made in other studies. This more realistic deduced overall heat loss coefficient will be compared to the one obtained by the Klein empirical formula, [12]. SYSTEM THEORY Under steady-state conditions, the collector delivers a useful energy rate (Qu) equal to the rate of radiation energy absorbed by the collector minus the energy rate transmitted to its Equations (2.4) and (2.5) will let us draw the efficiency curve versus from experimental surroundings. ..... .................................. ...................... ....... (2.1) After some reshuffling, it becomes, [17]: ..... .......................................... ........ (2.2) with FR, the heat removal factor that is determined by: .......... ..... . ..................... (2.3) The collector efficiency is given by: ..... ..... .. ....... .................. ... . (2.4) Which is equivalent to: .......... .......... (2.5) ....... data. Transforming equation (2.2), we obtain: ............................................. (2.6) From equation (2.6), we obtain the temperature difference between the absorber temperature and the fluid inlet temperature: ............ .. .1....... (2.7) ...... Taking into consideration equation (2.4), we obtain: ...... ... (2.8) .............. From equations (2.4) and (2.7), we obtain: ...... ................ .............(2.9) . ..... Combining equation (2.9) and equation (2.4), we obtain:........................1........... (2.10) in which .......... (2.11) ..... From equation (2.10), we derive FR ...... ...... (2.12) ..... We observe that FR can be calculated based on physical parameters and operating parameters, including the absorber temperature, environment temperature and inlet fluid temperature. The collector efficiency coefficient is given by equation (2.4) and equation (2.5). ....... Considering the experimental data, we plot the efficiency coefficient versus. From the slope ....... ............ of the experimental curve .... versus we determine the overall experimental loss coefficient knowing that FR can be determined by equation (2.12). For comparison, we calculate the heat loss coefficient using the widely used empirical equation given by Klein, [12]: .................. ...... .. . . (2.13) .............................. ..... .. ............. with ...1.0.04h..0.0005h....1.0.091...... (2.14) ..365.9.1.0.00883.....0.000129...... (2.15) h......... (2.16) Due to a low bottom casing temperature, the radiative heat exchange with the surrounding can be neglected. Only the convection term is considered, [17]; thus, the heat transfer coefficient for the heat transmitted from the back to the surrounding is given by: ...... . (2.17) ..... . ..,... Similarly, the heat transfer coefficient for the heat transferred from the collector edges to the environment, still assuming that radiation heat transfer is negligible, is given by: ...... . (2.18) ..... . ..,... The overall heat loss coefficient for the collector is the sum of the heat loss coefficient for the top, the heat loss coefficient for the bottom, and the heat loss coefficient for the edges. ........................ (2.19) Regarding the flat plate solar water collector with air absorber cooling, we replace the previous expression by the sum of heat rate for both fluids, and Ti takes the smallest inlet temperature from both fluids. EXPERIMENTAL SET-UP To determine the collector characteristics experimentally, primarily the overall heat loss coefficient, an experimental bench was set up (Fig.1). In the first configuration, a cross-section of the solar air collector shows a glass cover, an air gap between the glass cover and the absorber plate, the absorber plate fins, the air channel under the absorber, back and edge insulation and the casing. In the second configuration, a cross-section of the solar water collector shows a glass cover, an air gap between the glass cover and the absorber plate, a grooved plate to which riser tubes are embedded, back and edge insulation and the casing. The third configuration is identical to the solar water collector to which an air channel is added between a finned absorber and the back insulation. The used cover glass has low iron content in the three cases and is sealed to the casing with windshield silicone sealer. The absorber plate of 1 mm thickness is made of aluminium in the three cases. The absorber coating is matt black painting. In the second and third configurations, the embedded tube risers are made of copper. The back and the edges are insulated using rigid polyurethane panels. The benches are instrumented with thermocouples type K to measure temperatures (cover, absorbing plate, inlet and outlet fluids, collector back, ambient, etc.). Air speed was measured using a hot wire CFM anemometer with 0.01m/s resolution. A Kipp and Zonen pyranometer mounted at the collector orientation is used to measure the global incident solar radiation. During the experiment, recorded parameters come under the required environment conditions (Table 1). Physical properties of the materials were obtained from their documentation. All geometrical dimensions are given in the figures and Table 2. Airflow is generated with a centrifugal exhaust fan. In the case of the water collector, a water pump is used to perform forced circulation in the solar collector circuit. A flow meter is inserted in the water circuit in order to record the flow rate. Figure 2 shows the solar collectors constructed and used for the experiment. Table 1: Required environmental conditions (ASHRAE 93) Variable Absolute limits . Total solar irradiance normal to sun (W/m2) . Diffuse fraction (%) . Wind speed (m/s) . Incidence angle modifier . 790 (minimum) . 20 (maximum) . 2.28 d): 18.5 GBq per year, . For noble gases: the annual limit is calculated from the dose limit (50 µSv per year) on the border of the restricted protective zone according to adopted models, . For 3H and 14C there were no explicit limits for released activities. These limits were in force until 2007, when the limits for liquid effluents were modified to allow for a longer (18 months, previously 12 months) fuel cycle in NEK. The annual limit for the 3H release was increased, while the limit for other liquid releases was decreased. The new (and still valid) limits for liquid effluents are: . For all radionuclides except 3H, 14C and dissolved gases: 100 GBq per year, 40 GBq in a calendar quarter, . For 3H: 45 TBq per year, . For 14C: no explicit limit. Since 2007, the limits for gaseous releases and liquid effluents are a part of the Radiological Effluent Technical Specifications (RETS), a document related to radiological safety, monitoring, effluents, and reporting requirements. This document must be approved by a regulatory body for a valid operating licence and must be continuously maintained to reflect possible modifications and compliance with legal requirements. We can see that “the history” of dose and radioactive release limits for NEK is quite short and, except for the adjustment related to the introduction of 18-month fuel cycle, there has been no change of limits from the beginning of NEK operation. The reason is not in the indifference and passivity of the authorities but the professional knowledge of involved experts and adoption of internationally approved approach during the initial licensing process and initial operation of the plant. The author aims to review available data on emissions from NEK and compare it with available data on doses of the population in the NEK neighbourhood. We made such a comparison ten years ago, [1], and we would like to update data and verify conclusions from the previous paper. MONITORING The NEK must provide proof that it is complying with the imposed restrictions. Therefore, it is obliged to measure the emissions of radionuclides and report them to the regulatory body, and also to provide an independent assessment of effects to the environment and population. Monitoring, therefore, consists of the measurements of the plant releases (a measurement of emissions), sampling, and measurements in the plant neighbourhood (a measurement of immissions in the environment), evaluation of these measurements, as well as total dose assessment for members of the public based on collected data and/or computational models. Objectiveness and validity of results have been ascertained with the involvement of independent and authorised organizations (experts) in the monitoring implementation, evaluation of data and dose assessment. Requirements for monitoring in environment and emissions reporting have been in force from the very beginning of the operation of NEK. The programme of sampling and measurements, which has been verified and approved annually, has been based on a generic (and extensive) programme in the relevant rules from 1986 and the rules from 2007. The main difference between these rules is not in the content of the monitoring programme itself but the introduction of the additional requirements related to the quality and reliability of measurements. The rules were updated in 2018 when Council Directive 2013/59/Euratom was adopted in Slovenian legislation. Monitoring in the environment covers a 12 km circle around the plant and extends 30 km downstream of the River Sava. It includes: -external dose measurements with passive and active detectors, -sampling and measurement of radioactivity in the air (aerosols and iodine), -sampling and measurement of radioactivity in the river Sava (water, sediments, fish), -sampling and measurement of soil, -sampling and measurement of drinking water (wells and water supplies), -sampling and measurement of atmospheric precipitations and deposits, -sampling and measurement of food and milk (locally produced). In Figure 1, sampling locations for monitoring are presented (from [2]). Comparison of this figure with the figure from the report, [3], from 1991 reveals that the number of sampling points has not changed substantially from the late 1980s, and that sampling positions are not changed for the majority of "old" locations. In addition to sampling, monitoring also includes external dose measurements with passive TL detectors and continuous active detectors (energy compensated GM tubes) for external dose rate measurements. While the active detectors are a part of the early warning system in Slovenia and were not operational until the late 1990s, the network of TL detectors around NEK (57 locations, up to 10 km from the plant) has just slightly been changed from the 1990s. In support of monitoring and emergency preparedness, one meteorological station is situated inside the plant and three in the vicinity of the plant. The data from these stations serve for computational modelling required in dose assessment. All data from monitoring are collected and evaluated through authorised organizations on a regular basis. Evaluated data from environment monitoring and data on NEK emissions are used for the dose assessment of members of the public in the neighbourhood of NEK. Evaluated monitoring data, as well as results of dose assessment, are published on an annual basis as a report, which is also publicly available through the NEK internet site. Data from continuous detectors in the early warning system from all parts of Slovenia are publicly available on the internet, [4]. Figure 1: Current map of NEK neighbourhood with the distribution of places where different samples are regularly taken as a part of NEK monitoring, [2]. The position of NEK is marked with an arrow. OVERVIEW OF EMISSIONS FROM NEK As stated, for the purpose of monitoring, NEK provides data on released activities in all liquid effluents and gaseous releases from the beginning of operation (first criticality). All releases are continuously measured and data are collected. Released activity depends on many factors, but the most important factors are fuel quality, leaks from the primary system, produced power, trips and power excursions, the efficiency of liquid and gaseous waste processing, the chemistry of the primary system, and also an implementation of outages, and remediation of possible equipment failures. 4.1 Liquid effluents The history of annual activity releases in liquid effluents of NEK from 1983 to 2019 (data from [5]) is presented in Figure 2. We can see that the annual releases of fission and activation products (Figure 2a) have been steadily decreasing since the first years of operation. However, even the highest value (13.4 GBq in 1985) was a mere 6.7% of the annual limit (200 GBq per year until 2007). The values in 2008 and 2009 (the last one not presented in the figure) were under 100 MBq (less than 0.1% of the valid annual limit 100 GBq per year). The main contributors to annual releases were (and still are) activation products such as 58Co and 60Co (Figure 2c), and fission products 134Cs and 137Cs (Figure 2d). The released activity of 131I (Figure 2b) was always at least ten times lower and, in the last two decades, appears in the liquid effluents only occasionally. While the annual released activity of almost all radionuclides in liquid effluents has decreased substantially since the beginning of the operation, this does not apply to 3H (Figure 2e). The reason is in a direct connection between 3H production and energy production in NEK. The released activity was always close (40% to 80%) to the valid limit (20 TBq until 2007). In connection with the introduction of the 18-month fuel cycle, the Slovenian Nuclear Safety Administration has changed the limit to 45 TBq. 4.2 Gaseous releases The history of annually released activities in gaseous releases of NEK from 1983 to 2019 (data from [5]) is presented in Figure 3. The most important radionuclides are noble gases (133Xe, 133mXe, 41Ar), which are limited indirectly through effective dose on the NEK fence, and 3H (no limit). The highest releases of noble gases (Figure 3a) were in the mid-1990s (up to 25 TBq), while the current values are around 1 TBq. Similarly as in liquid effluents, an annual release of 3H in gaseous releases does not show any decrease (Figure 3b). After the introduction of an 18-month fuel cycle, the annually released activity is around 2.5 TBq. History of releases of 14C (Figure 3c), in contrast, shows an initial decrease from the first value of 400 GBq (1985) in the 1980s but has not substantially changed in the last two decades. The annual release is around 100 GBq. Annual releases of 131I (not presented in Figure 3) have been under 1 MBq since 2001. In the 1990s, these values were much higher, up to 2.7 GBq in 1996, which is 15% of the limiting value (18.5 GBq per year). In 1996, to our knowledge, was also the highest release activity of aerosols (20 MBq, 0.1% of limiting value). Annual activity releases of aerosols in the last decade are about MBq or less. Figure 2a: Fission and activation products. Figure 2b: Iodine (131I equivalent). The release The annual limit (200 GBq) was changed to is limited by the total limit for fission and 100 GBq 2007. activation products. Figure 2e: 3H (Tritium). The annual limit (20 GBq) was changed to 45 GBq in 2007. Figure 2: Annual activity releases in liquid effluents of NEK for 1983-2019 (from [5]). Annual releases of 14C are presented in Figure 3c. Releases were much higher in the 1980s than in the 1990s and later, when annual gaseous releases were 0.1 TBq per year or less. An annual limit has not been established for this radionuclide (the same also applies to other nuclear power plants) and the impact on population doses was not properly evaluated until 2002. Figure 3a: Noble Gases. The release is limited Figure 3b: 3H (Tritium). The annual limit has indirectly through a limited effective dose not been established. (50 µSv per year) on the border of the restricted area. Figure 3c: 14C. The annual limit has not been established. Figure 3: Annual activity releases in gaseous releases of NEK from 1983 to 2019 (from [5]). COMPARISON OF ASSESSED DOSES Data from the monitoring of emissions from NEK and immissions of radionuclides into the environment serve as a starting point for the evaluation of doses to the general population in the neighbourhood of NEK. This assessment is done by authorised organisations and finally submitted to the regulatory body through NEK. We will not discuss in detail methods or methodology of dose assessment in the past but merely present the data compiled from the available report on the evaluation of the radiological monitoring of the Krško Nuclear Power Plant in the 1980s, [3], and available national annual reports on the radiation and nuclear safety in the Republic of Slovenia, [6]. 5.1 Doses from liquid releases In Figure 4, doses from NEK liquid releases are presented with total doses from all releases (including gaseous). We can see that until 1994 almost all dose from NEK was associated with liquid releases. It was related to local fishermen and their families, and the most important radionuclide for this exposure pathway is 137Cs. However, we must point out that the total dose is not always a realistic term. It is not always possible to simply add doses from different exposure pathways since the complete dose assessment includes dose estimates for different and distinct critical groups of people. Comparison with data on Figure 2 can therefore hardly support the observable increase of annual doses in the first half of the 1990s in Figure 4. Doses from liquid releases were practically constant through the second half of the 1990s and diminished for almost an order of magnitude in 2002, which also applies to the total dose from NEK. Both changes were the consequence of the changed methodology for dose assessment. If we consider the last decade, the contribution of 14C to doses from liquid releases and to total doses was prevailing from 2013 to 2016 but became minimal in 2017. Figure 4: Annual doses of the population in the NEK neighbourhood due to liquid releases ([3], [6]). 5.2 Doses from gaseous releases Available data on annual doses from gaseous releases from NEK are presented in Figure 5. After 1994 (and until 2001), the most important exposure pathway was external exposure from the cloud of noble gases. Data in Figure 3a reveals that the release of noble gases was high in 1993 to 1996, but not in 1997 and afterward. However, assessed doses from noble gases were high until 2001. After 2001, the most important pathway was the ingestion of radionuclides from gaseous releases transferred to vegetation. The most important radionuclide in this respect is 14C, which becomes a part of the food chain. 14C was not included in the evaluation before 2002, and its contribution was not regularly measured until the last decade. Contribution to the dose in 2002 (and in the following years) was estimated based on analogous results with other nuclear power plants. This was also the consequence of the methodology change mentioned in the discussion of the doses from liquid releases. In this period, the inhalation dose from 3H was comparable to the estimated dose from 14C. In the previous decade, the doses from gaseous releases were related mostly to 14C. These doses were assessed from measurements in the environment and could be considered realistic values for the NEK neighbourhood. Figure 5: Annual doses of the population in the NEK neighbourhood due to gaseous releases ([3], [6]). 5.3 Importance of pathways In Figure 6, doses from liquid and gaseous releases are presented together with estimated total doses from NEK releases. The scale in the figure is linear to enable easier comparison of release pathway importance through the period from 1981 to 2019. In Figure 6, we can see that until 1994, total annual doses were approximately 2 to 4 µSv per year. The most important contribution was the ingestion dose due to liquid effluents. From 1995 to 2001, total doses were from 5 to 15 µSv per year, but the most important contribution comes from external exposure to gaseous releases (up to 11 µSv per year). This was the external exposure to radiation from a cloud of noble gases released from the plant. The assessed ingestion dose due to liquid releases was still high (2 to 6 µSv per year). From 2002 to 2019, total doses are (less than) 1 to (less than) 2 µSv per year. The most important contribution is the ingestion dose due to the release of 14C. The estimate considers the transfer of 14C from air to local vegetation and food chain and 14C in liquid releases which is consumed by local fishermen. All other pathways contribute much lower doses. We can also see that estimated doses have never come close to the authorised limit. Even the highest estimated dose (from 1995) represents only 30% of the authorised limit. Figure 6: Annual doses of the population in the NEK neighbourhood due to all releases ([3], [6]). The scale is linear to enable easier comparison of release pathways importance. 5.4 Comparison with natural background and global contamination Measurements of radioactivity in the neighbourhood of NEK provide us also with the data that enable us to assess natural radiation background. Similar data are available also for Ljubljana and other parts of Slovenia. Figure 7 presents annual external exposures to radiation measured with the TLDs in the neighbourhood of NEK, average external exposures in Slovenia (also measured with TLDs) and estimated annual external doses due to contamination from nuclear tests and Chernobyl releases in the NEK neighbourhood. The highest contribution of contamination to external dose was in 1986 (after Chernobyl) when the annual external dose from contamination was estimated to 325 µSv. We can see that the annual external dose, as measured with the TLDs in the NEK neighbourhood, was regularly slightly below the Slovenian average. We can also see that the contribution of external radiation to dose was significant only in the 1980s. In Figure 8, the total annual doses to the population in the NEK neighbourhood from NEK operation are compared to total doses from natural sources (natural background) and total doses due to contamination from nuclear tests and Chernobyl releases. We can see that measured natural background (2440 to 2530 µSv per year, which is practically the same as measured in Ljubljana) was always the main contributor to the total dose of population. In 1986, when the Chernobyl accident happened, the contribution of total contamination was 570 µSv (in Ljubljana even 720 µSv!), while is the estimated current value is about 10 µSv per year. The highest estimated annual dose from NEK (15.4 µSv in 1995) was 31% of the administrative dose limit and only 0.6% of natural background, while the (overestimated) annual doses from recent years (less than 0.15 µSv) present only 0.3% of administrative dose limit and 0.006% of natural background. Figure 7: Comparison of average annual external doses measured in the neighbourhood of NEK (57 TLDs), average external doses measured in Slovenia to annual external doses contamination estimated from measurements of soil contamination in NEK neighbourhood ([3], [6]). Figure 8: Comparison of annual doses of the population in the neighbourhood of NEK from exposure to natural sources, general contamination from nuclear tests and Chernobyl releases, and nuclear power plant. CONCLUSIONS We have seen that initial limitations imposed on the operation of NEK remain valid and, except in the case of liquid discharge of 3H, there was no need to modify them. NEK successfully operates within these limitations with releases, which are well below the limiting values. The scope and extent of monitoring have not changed substantially from the mid-1980s; therefore, old data remain compatible with recent measurements. We can recognise different approaches to dose evaluation in certain periods of monitoring history. Therefore, doses from different periods are not always directly comparable. According to the available data, except in the late 1990s, estimated doses to the population were always only a few percent points of the authorised dose limit, which is 50 µSv per year. References [1] M. Koželj: A brief history of Krško NPP radiation impact on environment, Proceedings of the International Conference Nuclear Energy for New Europe 2010, Portorož, Slovenia, September 6-9, 2010 [2] Nuklearna Elektrarna Krško: Meritve emisij in radioaktivnosti v okolju, Zemljevid merilnih mest [online]. Available: https://www.nek.si/sl/okolje/meritve-emisij-in­radioaktivnosti-v-okolju (10 February 2021) [3] B. Breznik, A. Kovac: Porocilo desetletnega obdobja radiološkega nadzora v okolici Nuklearne elektrarne Krško, Nuklearna elektrarna Krško, Krško, 1991 [4] Uprava Republike Slovenije za jedrsko varnost: Mreža zgodnjega opozarjanja, Trenutne vrednosti sevanja gama v okolju [online], Radioaktivnost v okolju. Available: https://radioaktivnost.si/rvo_public/RVO/Map (15. February 2021) [5] Uprava Republike Slovenije za jedrsko varnost: Razširjeno porocilo o varstvu pred ionizirajocimi sevanji in jedrski varnosti v Republiki Sloveniji leta 2019, URSJV/DP­216/2020, Ed. Benja Režonja, Uprava RS za jedrsko varnost, 1920 [6] Uprava Republike Slovenije za jedrsko varnost: Annual Reports on Radiation Protection and Nuclear Safety in Slovenia [online], 1992-2020. Available: https://www.gov.si/drzavni-organi/organi-v-sestavi/uprava-za-varstvo-pred-sevanji/o­upravi-republike-slovenije-za-varstvo-pred-sevanji/letna-porocila-o-varstvu-pred­ionizirajocimi-sevanji-in-jedrski-varnosti-v-republiki-sloveniji/ (10 February 2021) 40 JET JET Volume 14 (2021) p.p. 41-55 Issue 1, April 2021 Type of article 1.01 www.fe.um.si/en/jet.html COMPARISON OF CAVITATION MODELS FOR THE PREDICTION OF CAVITATION AROUND A HYDROFOIL PRIMERJAVA KAVITACIJSKIH MODELOV ZA NUMERICNO NAPOVED KAVITACIJE NA HIDRODINAMICNEM PROFILU Marko PezdevšekR, Ignacijo Biluš, Gorazd Hren Keywords: Hydrofoil, cavitation, Ansys CFX Abstract In this paper, four different cavitation models were compared for predicting cavitation around a hy­drofoil. A blocked structured mesh was created in ICEM CFD. Steady-state 2D simulations were per­formed in Ansys CFX. For all cases, the SST turbulence model with Reboud’s correction was used. For Zwart and Schnerr cavitation models, the recommended values were used for the empirical co­efficients. For the full cavitation model and Kunz cavitation model, values for the empirical coeffici­ents were determined as the recommended values did not provide satisfactory results. For the full cavitation model, the effect of non-condensable gases was neglected. For all the above-mentioned cavitation models, the pressure coefficient distribution was compared to experimental results from the literature. Povzetek V prispevku je narejena primerjava med štirimi kavitacijskimi modeli pri numericni napovedi kavita­cije na hidrodinamicnem profilu. V ICEM CFD je bila izdelana blokovna strukturirana mreža. V Ansys CFX so se izvedle 2D stacionarne simulacije. Za vse simulacije je bil uporabljen SST turbulentni model s korekcijo, ki jo je uvedel Reboud. Za kavitacijska modela Zwart in Schnerr smo uporabili privzete vrednosti empiricni koeficientov. Za full cavitation model in Kunzov kavitacijski model smo vrednosti R Corresponding author: Marko Pezdevšek, University of Maribor, Faculty of Energy Technology, Hocevarjev trg 1, 8270 Krško, E-mail address: marko.pezdevsek@um.si koeficientov dolocili sami, saj privzete vrednosti niso dale zadovoljivih rezultatov. Za vse štiri zgoraj omenjene kavitacijske modele smo primerjali porazdelitev tlacnega koeficienta z eksperimentalnimi rezultati iz literature. INTRODUCTION Cavitation is a phenomenon that occurs when a combination of low local static pressure and high velocities leads to pressures lower than the vapour pressure. Vapour structures occur in locations where the local pressure is below the vapour pressure. In some areas, cavitation can be beneficial, for example, in the medical field to remove kidney stones, but in engineering applications such as turbines, pump, and rudders, it is an undesirable effect. Cavitation may cause deterioration in performance, vibrations, and noise. Cavitation erosion occurs when the cavities collapse near the surface of a blade. Cavitation erosion is usually combined with other before mentioned unwanted cavitation effects. The development of cavitation in liquids can take different patterns. Typical types of cavitation have been classified based on their physical appearance. Some typical cavitation types are presented in Figure 1 and include bubble cavitation, sheet cavitation, vortex cavitation and cloud cavitation. Figure 1: Typical cavitation types, upper left travelling bubble cavitation, upper right leading edge sheet cavity, lower left vortex cavitation and lower right cloud cavitation, [1]. To reduce the cost of maintenance and improve the overall performance of a turbine, propellers, pumps or similar machinery, understanding and predicting cavitation and its effects is crucial. GEOMETRY AND MESH The hydrofoil geometry was obtained from [2]. As seen in Figure 2, the chord length of the hydrofoil is 152.4 mm, and the angle of attack is 1 °. The size of the domain also shown in Figure 2 is four chord lengths before, six chord lengths after and 2.5 chord length below and above the hydrofoil. Figure 2: Model dimensions. For the hydrofoil domain, a blocked structured mesh was created in ICEM CFD. The final mesh consisted of approximately 76,500 elements. The maximum dimensionless value y+ is below 1, as shown in Figure 3. Figure 3: Dimensionless y+ values on the hydrofoil surface. The upper image in Figure 4 shows the surface mesh of the model. The middle image shows a magnified cut-out section of the domain, which shows the mesh distribution around the hydrofoil. The bottom image in Figure 4 is a cut out magnified section of the middle image, where mesh distribution near the hydrofoil surface is visible. Figure 4: Surface mesh of the hydrofoil domain (upper image), cut-out magnified section of the upper image (middle section) and cut-out magnified section of the middle image (lower image). GOVERNING EQUATIONS, CAVITATION MODELS In CFX, the homogenous mixture flow is governed by the following set of equations, phases are considered incompressible and share the same velocity field U: Continuity equation: ............ .1. 1. (3.1) ..... ..... Where: ....– time-averaged mixture velocity [m/s], ..... – interphase mass transfer rate due to cavitation [kg/młs], ..... – vapour density [kg/mł], ..... – liquid density [kg/mł]. Momentum equation for the liquid vapour mixture: .............. ....................................................... (3.2) ........ .. Where: ....– density of the water-vapour mixture [kg/mł], ....– time averaged pressure [Pa], ....– dynamic viscosity of the water-vapour mixture [kg/m s], ..... – turbulent viscosity [kg/m s]. Volume fraction equation for the liquid phase. ......... (3.3) ............................... Where: ....– water volume fraction [/]. The water volume fraction and vapour volume fraction are defined as: ..liquidvolume..vapour volume; (3.4) total volume total volume The relation between the water and vapour fraction can be expressed as: ....1 (3.5) The water-vapour mixture density can be defined as:................1........... (3.6) The water-vapour mixture dynamic viscosity is defined as:................1........... (3.7) 3.1 Turbulence model Two-equation turbulence models are very widely used, as they offer good compromises between numerical effort and computational accuracy. In the models, the velocity and length scale are solved using separate transport equations. The k-e and k-. two-equation models use the gradient diffusion hypothesis to relate the Reynolds stresses to the mean velocity gradients and the turbulent viscosity. The turbulent viscosity is modelled as the product of a turbulent velocity and turbulent length scale. The turbulence velocity scale is computed from the turbulent kinetic energy, which is provided from the solution of its transport equation. The turbulent length scale is estimated from two properties of the turbulence field, usually the turbulent kinetic energy and its dissipation rate. The dissipation rate of the turbulent kinetic energy is provided from the solution of its transport equation. 3.1.1 The Shear Stress Transport (SST) Model The SST turbulence model was proposed by Menter, [3], and is a blend between the k-. model for the region near the surface and k-e model for the outer region. The model consists of a transformation of the k-e model to a k-. formulation. This is achieved by the use of a blending function F1. F1 is equal to one near the surface and decreases to a value of zero outside the boundary layer. [1] The turbulent kinetic energy k is defined by: .............. ................. ........................... .............................. ........................................ (3.8) Where: ..... – production rate turbulence, ....– turbulent kinetic energy [m2/s2], ...... – buoyancy production term. The specific dissipation rate . is obtained: .............. ................ ........ ........................... .................... .......... ....................................... (3.9) Where: ....– specific dissipation rate [s-1], ...... – buoyancy term. The model constants are: ......0.09, ..5/9, .....0.075, ......2, ......2. If we use ....., ..... and ..... to represent the terms in the k-e, k-. and SST model then the coefficients of the SST model are a linear combination of the corresponding coefficients of the underlying models [1]: ..................1............ (3.10) The turbulent viscosity is modified to account for the transport of the turbulent shear stress. The turbulent viscosity is defined as: ......... ...... (3.11) ....max..........,.......... Where: ....– strain rate magnitude [s-1], ..... – constant (0.31), ..... – second blending function. 3.1.2 Reboud’s correction Two-equation turbulence models were developed for single phase flows; they tend to overestimate the turbulent viscosity in the region of transition between vapour and liquid phase and damp the unsteadiness of the cavitating regime, [1]. Rebound, [4], proposed a modification of the k-e turbulence model by reducing the turbulent viscosity in order to take into account the suggested two-phase flow effects on the turbulent structures, [1]. The density in the turbulent viscosity equation is now replaced with a density function and is written as: ................ . ...... ........ ...... .......... (3.12) Where: ..... – mixture density [kg/m3], ....– constant (10). 3.2 Cavitation models The specific interphase mass transfer rate m . was modelled using an appropriate cavitation model. We assume that the specific mass transfer rate is positive if directed from vapour to liquid. 3.2.1 Zwart The Zwart cavitation model was developed by Zwart et al. [5]. The model is based on the multiphase flow equations, with mass transfer due to cavitation appearing as source-and-sink terms in the liquid and vapour continuity equations. The mass transfer rate is derived from a simplified Rayleigh-Plesset model, [5]. 3........1........... .......... ........ .2 if.......... . ..... 3..... ..... . (3.13) . 3......... .......... . .........2 if .......... ..... 3...... Where: ....... – nucleation site volume fraction [m], ..... – bubble radius [m], ..... – vapour pressure [Pa], ....... – evaporation coefficient [/], ........ – condensation coefficient [/]. The recommended values for the two coefficients are ........50 and .........0.01. The recommended values for the nucleation site volume fraction and bubble radius are ........5· 10.. and ......10.. . 3.2.2 Schnerr Schnerr and Sauer, [6], assumed that the vapour structure is filled with spherical bubbles, which are governed by the simplified Rayleigh Plesset equation. The mass transfer rate in the Schnerr and Sauer model is proportional to .....1....... Moreover, the function .........1...... has the interesting property that it approaches zero when .....0and .....1and reaches the maximum in between, [1]. .......... .......... ............ .....1......3.2 if.......... . ..... 3..... ..... . (3.14) . ................... ...................1..........3..23.... if .......... .. Where: ..... – bubble radius [m], ....... – evaporation coefficient [/], ........ – condensation coefficient [/]. ....... .... 3 (3.15) 1.....4......... Where: ....– bubble number density [/], The recommended values for the two coefficients are ....... .1 and .........0.2. The recommended values for the bubble number density is .....10.. . 3.2.3 Full cavitation model The full cavitation model (FCM) was developed by Singhal et al. [7]. The bubble dynamics equation is referred to as a “reduced bubble dynamics formulation” and is derived from the generalized Rayleigh-Plesset equation. It assumes that in most engineering situations, there are plenty of nuclei for the inception of cavitation. The primary focus is on the proper account of bubble growth and collapse, [7]. v.............. .....................2.1............. if.......... . 3..... ..... . (3.16) . v.............. . ..... .......2..... if .......... .... ....3..... . Where: ..... – vapour mass fraction [/], ..... – gas mass fraction [/], ....– surface tension [N/m], ..... – vapour pressure [Pa], ..... – empirical coefficient [m/s], ..... – empirical coefficient [m/s]. Singhal et al., [8], reported a numerical model using a probability density function approach for accounting for the effects of turbulent pressure fluctuations. The local values of the turbulence pressure fluctuations were defined as: ......... .0.39........ (3.17) In [8], computations of time-averaged phase-change rates by the integration of instantaneous rates in conjunction with the assumed probability density function for pressure variation with time were presented. By raising the phase change threshold pressure value, Singhal, [7], proposed a simplified term: ...................2 . (3.18) where: ....... – saturation pressure [Pa], The recommended values for the two empirical calibration coefficients are ......0.02and ...... 0.01. 3.2.4 Kunz The Kunz cavitation model is a heuristic model based on work by Merkle et al., [9]. The source term is subdivided into a term related to vaporization and a term related to condensation. The transformation of liquid to vapour is calculated as proportional to the amount by which the pressure is below the vapour pressure. For the transformation of vapour to liquid, a simplified form of the Ginzburg-Landau potential is employed, [10]. . .......................... . ..... ...... .................min.0,..... (3.19) . . ..12................. where: ........ – empirical coefficient [/], ........ – empirical coefficient [/], ..... – free stream velocity [m/s], ..... – mean flow time scale [s]. The mean flow time scale is defined as: ........../..... (3.20) where: ....– characteristic length scale [m]. The recommended values for the two empirical coefficients are .........100and ........ .100. BOUNDARY CONDITIONS AND SIMULATION SETTINGS For this problem, the following boundary conditions were applied (shown in Figure 5): -The left surface was specified as an inlet, where a normal velocity of 16.91 m/s was defined. -The right surface was defined as an outlet, where a static pressure of 51,957 Pa was defined. -The top and bottom surfaces were defined as a free-slip wall. -The symmetry boundary condition was applied for the side surfaces. -The hydrofoil surface was defined as a no-slip wall. Figure 5: Boundary conditions. The cavitation number is a non-dimensional parameter for cavitating flow and is defined as: .. ............. (4.1) 0.5............. The cavitation number for this case is calculated to 0.34. For all four selected cavitation models, steady-state 2D simulations were performed in Ansys CFX. For all cases, the average RMS residuals were set at 10^-6, and the SST turbulence model with Reboud’s correction was used. For Zwart and Schnerr cavitation models, the recommended values were used for the empirical coefficients. For the FCM cavitation model, we have neglected the effect of non-condensable gases. The two empirical coefficients were set at ......1and ...... 1. For the Kunz cavitation model, the two empirical coefficients were set at ........ .65,000and .........1800. The empirical coefficients used for this study are presented in Table 1. Table 1: Empirical coefficient values for cavitation models used in this study. Cavitation model Coefficient values Zwart ....... .50, ........ .0.01 Schnerr ....... .1, ........ .0.2 FCM ..... .1, ..... .1 Kunz ........ .65,000, ........ .1800 RESULTS Figure 6 shows the distribution of the pressure coefficient along the hydrofoil surface for all four selected cavitation models. The pressure coefficient is compared to the experimental results from the literature, [2]. In general, all cavitation models show good agreement with experimental results at x/c values below 0.7. Above 0.7, the FCM and Zwart cavitation models show the best agreement, while the Kunz cavitation model deviates the most. Cp [/] -0,4 -0,3 -0,2 -0,1 0,0 0,1 0,2 0,3 x/c [/] Figure 6: Pressure coefficient distribution. The cavitation for the Schnerr model starts at x/c=0.2; for the FCM and Kunz model, it starts at approximately x/c=0.3. For the Zwart model, the cavitation starts at x/c=0.4. Vapour volume fraction distribution for all cavitation models is seen in Figure 7. 1,0 0,9 Volume vapour fraction [/] 0,8 0,7 0,6 0,5 0,4 0,3 0,2 0,1 0,0 x/c [/] Figure 7: Vapour volume fraction distribution. Figure 8 shows the vapour volume fraction for all cavitation models. The Zwart and FCM models predicted cavitation in a very similar manner, which is also evident from the pressure coefficient distribution. For the Schnerr cavitation model, the values for the vapour volume fraction are higher compared to the other three models. CONCLUSION Steady-state 2D simulations were performed for a hydrofoil with a chord length of 152.4 mm and an angle of attack of 1 °. A velocity of 16.91 m/s was defined, the cavitation number was calculated to 0.34. For all cases, the SST turbulence model with Reboud’s correction was used. For Zwart and Schnerr cavitation models, the recommended values were used for the empirical coefficients; the results for the pressure coefficient for both models show good agreement with experimental results. For the FCM cavitation model, we have neglected the effect of non-condensable gases. In this study, the two empirical coefficients were set at ......1and ......1. With the set coefficients, the results for the pressure coefficient were in good agreement with the experimental results. For the Kunz cavitation model, the two empirical coefficients were set at ........ .65,000 and ........ .1800. With the set coefficients, the results for the pressure coefficient were in good agreement with the experimental results but compared to other models; the Kunz cavitation model deviates the most. From the results given, it seems that the selected cavitation models in this study can offer similar levels of accuracy, although we should note that with the FCM and Kunz model, the recommended values for the empirical coefficients did not provide satisfactory results. References [1] D. Ziru LI: Assessment of Cavitation Erosion with a Multiphase Reynolds-Averaged Navier-Stokes Method, Wuhan University of Technology, Wuhan, P.R. China, 2012. [2] ANSYS, Inc.: ANSYS CFX tutorials, 2016 [3] F. R. Menter: Two-Equation Eddy-Viscosity Turbulence Models for Engineering Applications, AIAA Journal, Volume 32, No. 8, August 1994. [4] J. L. Reboud, B. Stutz, O. Coutier-Delgosha: Two phase flow structure of cavitation: experiment and modeling of unsteady effects, Proceedings of the 3rd International Symposium on Cavitation, Grenoble, France, 1998. [5] P. J. Zwart, A. G. Gerber, T. Belamri: A Two-Phase Flow Model for Predicting Cavitation Dynamics, ICMF 2004 International Conference on Multiphase Flow, Yokohama, Japan, 2004 [6] G. H. Schnerr: Physical and Numerical Modeling of Unsteady Cavitation Dynamics, ICMF­2001, 4th International Conference on Multiphase Flow, New Orleans, USA, 2001 [7] A. K. Singhal, H. Li, Y. Jiang: Mathematical Basis and Validation of the Full Cavitation Model, Journal of Fluids Engineering, 2002 [8] A. K. Singhal, N. Vaidya, A. D. Leonard: Multi-Dimensional Simulation of Cavitating Flows Using a PDF Model for Phase Change, ASME FED Meeting, Vancouver, Canada, 1997. [9] C. L. Merkle, J. Feng, P. Buelow: Computational modelling of the dynamics of sheet cavitation, Proceedings of 3rd International Symposium on Cavitation, Grenoble, France, 1998. [10] R. F. Kunz, D. A. Boger, D. R. Stinebring, T. S. Chyczewski, J. W. Lindau, H. J. Gibeling, S. Venkateswaran, T. R. Govindan: Preconditioned Navier-Stokes Method for Two-Phase Flows with Application to Cavitation Prediction, Computers & Fluids 29, 2000. Nomenclature (Symbols) (Symbol meaning) .... time-averaged mixture velocity .... interphase mass transfer rate due to cavitation ........ vapour density ........ liquid density .... density of the water-vapour mixture .... time-averaged pressure .... dynamic viscosity of the water-vapour mixture ........ turbulent viscosity .... water volume fraction ........ production rate turbulence .... turbulent kinetic energy ............ buoyancy production term .... specific dissipation rate ............ buoyancy term .... strain rate magnitude ........ second blending function ........ mixture density nucleation site volume fraction ................ ........ bubble radius ........ vapour pressure ................ evaporation coefficient condensation coefficient .................... ........ bubble radius .... bubble number density ........ vapour mass fraction ........ gas mass fraction .... surface tension ........ vapour pressure ........ empirical coefficient ........ empirical coefficient ' .................... local values of the turbulence pressure fluctuations ................ saturation pressure .................... empirical coefficient .................... empirical coefficient ..... free stream velocity ....8 mean flow time scale .... characteristic length scale .... cavitation number 56 JET JET Volume 14 (2021) p.p. 57-67 Issue 1, April 2021 Type of article 1.01 www.fe.um.si/en/jet.html MODELLING OF MAGNETIC REGENERATOR AND HEAT TRANSFER AGENT IN MICROCHANNELS MODELIRANJE MAGNETNEGA REGENERATORJA IN SREDSTVA ZA PRENOS TOPLOTE V MIKROKANALIH Dorin Botoc1,R,Bianca Eliza Oneata2, Ionut Rusu1, Alexandru Salceanu1 , Jurij Avsec3 Keywords: magnetic refrigeration, active magnetic regenerator, magnetocaloric material Abstract In this article, a brief introduction of conventional refrigeration is given, followed by a description and history of magnetic refrigeration. The active magnetic regenerator comprises 12 parallel plates of ma-gnetocaloric material (gadolinium (Gd)), through which circulates the heat transfer fluid (water, in this case). At both ends of the regenerator are the heat exchangers. The hot heat exchanger (HHEX) and the cold heat exchanger (CHEX) connects the fluid to the heat sources. The principle of operation of a magnetic refrigeration installation is based on exploiting the magnetocaloric effect from the materials that possess these properties (in this case, Gd). R Corresponding author: Dorin Botoc, Faculty of Electrical Engineering, Energetics and Applied Informatics, Gheorghe Asa-chi Technical University of Iasi, Romania, E -mail address: HYPERLINK "mailto:dorinbotoc@yahoo.com" dorinbotoc@ya-hoo.com 1Faculty of Electrical Engineering,Energetics and Applied Informatics,Gheorghe Asachi Technical University of Iasi, Romania 2 University of Strasbourg, Foreign Applied Languages: English - French and Spanish, France 3 University of Maribor, Faculty of Energy Technology, Laboratory for Thermomechanics, Applied Thermal Energy Technolo­gies and Nanotechnologies, Krško, Slovenia Povzetek V tem clanku je uvodoma opisano kompresorskego hlajenje, ki mu sledi opis in zgodovina magnetnega hlajenja. Aktivni magnetni regenerator obsega 12 vzporednih plošc magnetokaloricnega materiala (gadolinij (Gd)), skozi katerega kroži tekocina za prenos toplote (v tem primeru voda). Na obeh koncih regeneratorja sta izmenjevalnika toplote. Vroci izmenjevalnik toplote (HHEX) in hladni izmenjevalnik toplote (CHEX) povezujeta tekocino z viri toplote. Nacelo delovanja magnetne hladilne naprave temelji na izkorišcanju magnetokaloricnega ucinka materialov, ki imajo te lastnosti (v tem primeru Gd). BACKGROUND AND INTRODUCTION The technology for producing artificial cold based on vapour compression was introduced more than 120 years ago, with small improvements since then, [1]. However, it has reached a technological level at which improvements and performance increases are unlikely. In 2002, they represented 25% of the total electricity consumption in residential areas and 15% of commercial electricity consumption in the US. A difficult problem to solve and quite serious with this type of technology is the refrigerant used; in most cases, it was hydrofluorocarbons (HFC), which when removed from use; this compressible gas is very harmful both for the environment as well as for the ozone layer, [2]. There is a great deal of research on the evolution and development of freons with a much lower pollution level and contribution to global warming. Previously used refrigerants, specifically chlorofluorocarbons (CFCs) and hydrochlorofluorocarbons (HCFCs), have had a devastating effect on the ozone layer, [3]. In this century, with the increased concern regarding global warming and increasing energy consumption, the development of a more energy-efficient refrigeration technology than steam compression is a priority, [4]. Meeting the growing demands of refrigeration and air conditioning while reducing greenhouse gases has led to further research and developments in magnetic refrigeration systems to produce artificial refrigeration. This motivation is the basis for attracting interest in many research centres in Europe, the United States, China, and Japan. Magnetic refrigeration is based on the magnetocaloric effect of materials with this property, which manifests itself as a change in temperature when subjected to a magnetic field, [5]. This effect is reversible: when the material is no longer under the action of the magnetic field, it returns to its original state. The magnitude of the temperature change depends, in most cases, on the strength of the magnetic field. To obtain maximum energy efficiency and minimize the operating costs of this installation, it is recommended that the magnetic field source use permanent magnets, [6]. As a working agent, magnetic refrigeration uses solid materials: silicon compounds, gadolinium, etc. These materials illustrate the unique and most strongly highlighted property: the magnetocaloric effect, which is manifested by the increase/decrease of the temperature of that material when it is magnetized/demagnetized. This process of obtaining low temperatures (i.e., lower than those of the environment) have thermal conditioning as the main purpose. BOUNDARY CONDITION The modelling and simulation of thermodynamic processes were performed in COMSOL Multiphysics software. The model was inspired by specialized publications and was built using the options offered by this software. The simulated model is assumed to have a mesh model. In Fig. 1, the AMR (active magnetic regenerator) model is presented. Figure 1: COMSOL Multiphysics active magnetic regenerator model used for the simulation Heat exchangers (HHEX and CHEX) are located at both ends of the AMR. The fluid channel length of the modelled AMR is 300 mm, and the width is 10 mm. The materials were also selected according to the simulations performed in the profile publications, [7]. The AMR uses magnetocaloric gadolinium (Gd) and the assembly was made in the form of parallel plates. Water was used as a heat transfer fluid, [8]. Figure 2: Mesh of the fluid channel modelled in COMSOL Multiphysics MODELLING AND SIMULATION This modelling and simulation are mainly based on the study of fluid channels inside the AMR under different magnetization conditions. The modelling was done in COMSOL Multiphysics, and the aim is to increase the efficiency of the regenerator with 12 parallel plates and to use it as a heat-transfer agent. Figure 3: 3D model of the fluid channel modelled in COMSOL Multiphysics without Gd plates The following are the geometric features of this model. The geometry is with parallel plates and microchannels through which the flow of the working agent occurs, [9]. Figure 4: 2D model microchannels for heat transfer agent, modelled in COMSOL Multiphysics with Gd plates The heat transfer fluid used in this AMR model is water, and COMSOL Multiphysics has already implemented Newtonian fluid in its library. To avoid the oxidation of gadolinium (Gd) plates over time, various compounds can be added to the working agent (water) provided that the physical and thermodynamic properties are not significantly altered. There are different anticorrosive compounds for these types of installations. Figure 5: 3D model of fluid flow along the length of the Gd plate. inside the canal When studying the fluid flow, the first characteristic criterion to be considered is the Reynolds number in order to determine the flow regime, [10]. ................. (3.1) ..... Where: .....,....,..... are, density, speed and dynamic viscosity of the heat transfer fluid, and ..... is the hydraulic diameter. The hydraulic diameter represents the characteristic length of the fluid flow and can be defined as: (3.2) ......4......... Where: ..... is the flow section and ...., wetted perimeter RESULTS AND DISCUSSION When the heat transfer agent flows in the X direction, having been in contact with the Gd plates, which have a different temperature, the so-called hydrodynamic boundary layer appears near the solid wall. The profile below is circular, so it is possible to define an average Reynolds number, [11-12]. It can be seen how the fluid temperature drops from a maximum in the centre of the microchannel to almost zero to the boundary layer of the channel wall. Figure 6: 2D model of the fluid channel temperature modelled in COMSOL Multiphysics Magnetization, specific heat and adiabatic variation of the temperature of the magnetocaloric material Gd varies depending on two essential parameters: temperature and the applied magnetic field, [13]. Figure 7: 3D model of the fluid channel temperature modelled in COMSOL Multiphysics with 12 parallel Gd plates In this way, a direct application of the adiabatic temperature variation on the parallel Gd plates is obtained. At each step, the temperature of the plates made of magnetocaloric material Gd. varies according to the equation below, [14]: ......................,...............,...................... (4.1) The internal energy for the plates of magnetocaloric material Gd can be written as follows:.....,......... (4.2) And in differential form: .................... (4.3) For the operation of the magnetic refrigeration system at constant pressure, the enthalpy can be defined as: .................. (4.4) And in differential form: .................... (4.5) Fig.9 shows the variation of the adiabatic temperature of the gadolinium for a variation of the magnetic field from 2T (black line) to the maximum intensity of the magnetic field of 10T (light blue line). We can see the evolution of the values of the adiabatic temperature variation of gadolinium as a function of temperature. Figure 8: Variation of the adiabatic temperature of the Gd plates made with the MatLab algorithm In Fig. 10, we can see the evolution of the entropy change of the Gd plates depending on the applied magnetic field. This variation of the magnetic field is made from 2T (black curve) to 10T (light blue curve). For all values of magnetic intensity, the evolution is slow at the beginning until the maximum value of 15 J/kgK, after which it is slowly decreasing. We observed that the contribution of electron entropy is negligible, and after Curie temperatures, the magnetocaloric material gadolinium becomes paramagnetic, [15]. Figure 9: Evolution of the entropy change of Gd plates. performed in the MatLab algorithm In Fig. 10 (below), one can see the 3D model simulated in COMSOL Multiphysics of the fluid channel in AMR with 12 gadolinium plates (Gd.) CONCLUSIONS AND PERSPECTIVES The first part of this paper presented the basic concepts of magnetic refrigeration and its definitions. In recent years, interest in this technology has grown considerably, due to the advantages we have shown. Then, we explained how this technology works based on magnetocaloric materials, in our case, gadolinium (Gd), which remains the reference material in this field due to its specific properties. In addition to the magnetocaloric material that is part of the AMR, the heat transfer fluid and the location of the parallel plates are quite important components. In Fig. 1 an AMR with 12 parallel gadolinium plates (Gd) was modelled. Fig. 2 shows the fluid channel mesh; 12 parallel gadolinium plates (Gd) were used. Fig. 10 shows the 3D model of the 12 gadolinium plates (Gd) in COMSOL Multiphysics, placed in parallel and water was used as a working agent. The advantages of this technology were shown with the possibility of further study of different types of AMRs: with gadolinium powder (Gd), with gadolinium spheres (Gd) and the accentuation of the modelling of heat transfer, magnetism, and convective transport of the working agent. References [1] Sergiu Lionte: Caractérisation, étude et modélisation du comportement thermomagnétique d’un dispositif de réfrigération magnétique ŕ matériaux non linéaires et point de Curie proche de la température ambiante, PhD Thesis, 2015 [2] Brown, G. V.:Magnetic heat pumping near room temperature, Journal of Alloys and Compounds, 47(8), 3673-3680, 1976 [3] Kitanovski, A. Egolf: Thermodynamics of magnetic refrigeration, International journal of refrigeration, 29(3-12), 3-21, 2005 [4] J Tusek, A Kitanovski, A Flisar et al:Active magnetic regenerator (AMR) experimental devices, în Fifth IIF-IIR international conference on magnetic refrigeration at room temperature, Thermag V, Grenoble, France, 2012 [5] Nellis, G., Engelbrecht, K., Klein, S., Zimm, C., and Russek, S.:Model and Evaluation Tools for Assessing Magnetocaloric Effect for Space Conditioning and Refrigeration Applications, 2004 [6] Tishin, A. M., Gschneidner, K. A., and Pecharsky, V. K.:Magnetocaloric effect and heat capacity in the phase-transition region, Physical Review B, 59(1), 503-511, 1999 [7] Tishin, A., and Spichkin, Y.:The Magnetocaloric Effect and its Applications., Institute of Physics Publishing, 2003 [8] Kirol LD, Mills JI:Numerical analysis of thermomagnetic generators, J Appl Phys 56:824–828, 1984 [9] Salomon D:Thermomagnetic mechanical heat engines, J Appl Phys 65:3687–3693, 1989 [10] Pecharsky VK, Gschneidner KA Jr:Effect of alloying on the giant magnetocaloric effect of Gd5(Si2Ge2), J Magn Magn Mater 167:179–184, 1997 [11] Sergiu Lionte et al.:Numerical analysis of a reciprocating active magnetic regenerator, Applied Thermal Engineering, ELSEVIER, 2014 [12] C Vasile, C Muller:Innovative design of a magnetocaloric system, Int J Refrig, nr. 29:1318-1326, 2006 [13] C Vasile, C Muller:A new system for a magnetocaloric refrigerator, în First IIF-IIR international conference on magnetic refriation at room temperature, Thermag I, Montreux, Switzerland, 2005 [14] KW Lipso, KK Nilesen, DV Christensen et al: Measuring the effect of demagnetization in stacks of gadolinium plate using the magnetocaloric effect, J Magn Magn Mater, nr. 323:3027-3032, 2011 [15] C Aprea, A Grego, A Maiorino et al:Initial experimental results from a rotary permanent magnet magnetic refrigerator,Int J Refrig, vol. 43, nr. 111-122, 2014 Nomenclature (Symbols) (Symbol meaning) AMR active magnetic regenerator Gd gadolinium material HHEX hot heat exchanger CHEX cold heat exchanger 68 JET Author instructions www.fe.um.si/en/jet.html MAIN TITLE OF THE PAPER SLOVENIAN TITLE Author1, Author 2, Corresponding authorR Keywords: (Up to 10 keywords) Abstract Abstract should be up to 500 words long, with no pictures, photos, equations, tables, only text. Povzetek (Abstract in Slovenian language) Submission of Manuscripts: All manuscripts must be submitted in English by e-mail to the editorial office at jet@um.si to ensure fast processing. Instructions for authors are also available online at http://www.fe.um.si/en/jet/author-instructions.html. Preparation of manuscripts: Manuscripts must be typed in English in prescribed journal form (MS Word editor). A MS Word template is available at the Journal Home page. A title page consists of the main title in the English and Slovenian language; the author(s) name(s) as well as the address, affiliation, E-mail address, telephone and fax numbers of author(s). Corre­sponding author must be indicated. Main title: should be centred and written with capital letters (ARIAL bold 18 pt), in first paragraph in English language, in second paragraph in Slovenian language. Key words: A list of 3 up to 6 key words is essential for indexing purposes. (CALIBRI 10pt) Abstract: Abstract should be up to 500 words long, with no pictures, photos, equations, tables, ­ text only. Povzetek: - Abstract in Slovenian language. Main text should be structured logically in chapters, sections and sub-sections. Type of letters is Calibri, 10pt, full justified. R Corresponding author: Title, Name and Surname, Organisation, Department, Address, Tel.: +XXX x xxx xxx, E-mail address: x.x@xxx.xx 1 Organisation, Department, Address 2 Organisation, Department, Address Units and abbreviations: Required are SI units. Abbreviations must be given in text when first men­tioned. Proofreading: The proof will be send by e-mail to the corresponding author in MS Word’s Track changes function. Corresponding author is required to make their proof corrections with accepting or rejecting the tracked changes in document and answer all open comments of proof reader. The corresponding author is responsible to introduce corrections of data in the paper. The Editors are not responsible for damage or loss of submitted text. Contributors are advised to keep copies of their texts, illustrations and all other materials. The statements, opinions and data contained in this publication are solely those of the individual authors and not of the publisher and the Editors. Neither the publisher nor the Editors can accept any legal responsibility for errors that could appear during the process. Copyright: Submissions of a publication article implies transfer of the copyright from the author(s) to the publisher upon acceptance of the paper. Accepted papers become the permanent property of “Journal of Energy Technology”. All articles published in this journal are protected by copyright, which covers the exclusive rights to reproduce and distribute the article as well as all translation rights. No material can be published without written permission of the publisher. Chapter examples: MAIN CHAPTER (Arial bold, 12pt, after paragraph 6pt space) 1.1 Section (Arial bold, 11pt, after paragraph 6pt space) 1.1.1 Sub-section (Arial bold, 10pt, after paragraph 6pt space) Example of Equation (lined 2 cm from left margin, equation number in normal brackets (section. equation number), lined right margin, paragraph space 6pt before in after line): Equation (1.1) Tables should have a legend that includes the title of the table at the top of the table. Each table should be cited in the text. Table legend example: Table 1: Name of the table (centred, on top of the table) Figures and images should be labelled sequentially numbered (Arabic numbers) and cited in the text – Fig.1 or Figure 1. The legend should be below the image, picture, photo or drawing. Figure legend example: Figure 1: Name of the figure (centred, on bottom of figure, photo, or drawing) References [1] N. Surname: Title, Journal Title, Vol., Iss., p.p., Year of Publication [2] N. Surname: Title, Publisher, Year of Publication [3] N. Surname: Title [online], Publisher or Journal Title, Vol., Iss., p.p., Year of Publication. Available: website (date accessed) Examples: [1] J. Usenik: Mathematical model of the power supply system control, Journal of Energy Technology, Vol. 2, Iss. 3, p.p. 29 – 46, 2009 [2] J. J. DiStefano, A.R. Stubberud, I. J. Williams: Theory and Problems of Feedback and Con­trol Systems, McGraw-Hill Book Company, 1987 [3] T. Žagar, L. Kegel: Preparation of National programme for SF and RW management taking into account the possible future evolution of ERDO [online], Journal of Energy Technol­ogy, Vol. 9, Iss. 1, p.p. 39 – 50, 2016. Available: http://www.fe.um.si/images/jet /Volume 9_Issue1/03-JET_marec_2016-PREPARATION_OF_NATIONAL.pdf (7. 10. 2016) Example of reference-1 citation: In text [1], text continue. Nomenclature (Symbols) (Symbol meaning) t time ISSN 1855-5748 9 771855 574008 JET I Journal of Energy Technology I Vol. 14, Issue 1, April 2021 I University of Maribor, Faculty of Energy Technology