„ Since 1955 . Strojniski vestnik Journal of Mechanical Engineering Strojniški vestnik - Journal of Mechanical Engineering (SV-JME) Aim and Scope The international journal publishes original and (mini)review articles covering the concepts of materials science, mechanics, kinematics, thermodynamics, energy and environment, mechatronics and robotics, fluid mechanics, tribology, cybernetics, industrial engineering and structural analysis. The journal follows new trends and progress proven practice in the mechanical engineering and also in the closely related sciences as are electrical, civil and process engineering, medicine, microbiology, ecology, agriculture, transport systems, aviation, and others, thus creating a unique forum for interdisciplinary or multidisciplinary dialogue. The international conferences selected papers are welcome for publishing as a special issue of SV-JME with invited co-editor(s). Editor in Chief Vincenc Butala University of Ljubljana, Faculty of Mechanical Engineering, Slovenia Technical Editor Pika Škraba University of Ljubljana, Faculty of Mechanical Engineering, Slovenia Founding Editor Bojan Kraut University of Ljubljana, Faculty of Mechanical Engineering, Slovenia Editorial Office University of Ljubljana, Faculty of Mechanical Engineering SV-JME, Aškerčeva 6, SI-1000 Ljubljana, Slovenia Phone: 386 (0)1 4771 137 Fax: 386 (0)1 2518 567 info@sv-jme.eu, http://www. sv-jme.eu Print: Papirografika, printed in 300 copies Founders and Publishers University of Ljubljana, Faculty of Mechanical Engineering, Slovenia University of Maribor, Faculty of Mechanical Engineering, Slovenia Association of Mechanical Engineers of Slovenia Chamber of Commerce and Industry of Slovenia, Metal Processing Industry Association President of Publishing Council Mitjan Kalin University of Ljubljana, Faculty of Mechanical Engineering, Slovenia Vice-President of Publishing Council Bojan Dolšak University of Maribor, Faculty of Mechanical Engineering, Slovenia Cover: Time-frequency maps of the modelled and the experimental pendant droplet vertical position in the process of continuous droplet generation with mass spring resonant detachment regime. Examples of 2D and 3D deposited droplets. Image courtesy: University of Ljubljana, Faculty of Mechanical Engineering, Laboratory of Synergetics ISSN 0039-2480, ISSN 2536-2948 (online) International Editorial Board Kamil Arslan, Karabuk University, Turkey Hafiz Muhammad Ali, University of Engineering and Technology, Pakistan Josep M. Bergada, Politechnical University of Catalonia, Spain Anton Bergant, Litostroj Power, Slovenia Miha Boltežar, University of Ljubljana, Slovenia Filippo Cianetti, University of Perugia, Italy Franci Čuš, University of Maribor, Slovenia Janez Diaci, University of Ljubljana, Slovenia Anselmo Eduardo Diniz, State University of Campinas, Brazil Jožef Duhovnik, University of Ljubljana, Slovenia Igor Emri, University of Ljubljana, Slovenia Imre Felde, Obuda University, Faculty of Informatics, Hungary Janez Grum, University of Ljubljana, Slovenia Imre Horvath, Delft University of Technology, The Netherlands Aleš Hribernik, University of Maribor, Slovenia Soichi Ibaraki, Kyoto University, Department of Micro Eng., Japan Julius Kaplunov, Brunel University, West London, UK Iyas Khader, Fraunhofer Institute for Mechanics of Materials, Germany Jernej Klemenc, University of Ljubljana, Slovenia Milan Kljajin, J.J. Strossmayer University of Osijek, Croatia Peter Krajnik, Chalmers University of Technology, Sweden Janez Kušar, University of Ljubljana, Slovenia Gorazd Lojen, University of Maribor, Slovenia Thomas Lubben, University of Bremen, Germany Jure Marn, University of Maribor, Slovenia George K. Nikas, KADMOS Engineering, UK Tomaž Pepelnjak, University of Ljubljana, Slovenia Vladimir Popovič, University of Belgrade, Serbia Franci Pušavec, University of Ljubljana, Slovenia Mohammad Reza Safaei, Florida International University, USA Marco Sortino, University of Udine, Italy Branko Vasič, University of Belgrade, Serbia Arkady Voloshin, Lehigh University, Bethlehem, USA General information Strojniški vestnik - Journal of Mechanical Engineering is published in 11 issues per year (July and August is a double issue). Institutional prices include print & online access: institutional subscription price and foreign subscription €100,00 (the price of a single issue is €10,00); general public subscription and student subscription €50,00 (the price of a single issue is €5,00). Prices are exclusive of tax. Delivery is included in the price. The recipient is responsible for paying any import duties or taxes. Legal title passes to the customer on dispatch by our distributor. Single issues from current and recent volumes are available at the current single-issue price. To order the journal, please complete the form on our website. For submissions, subscriptions and all other information please visit: http://www.sv-jme.eu. You can advertise on the inner and outer side of the back cover of the journal. The authors of the published papers are invited to send photos or pictures with short explanation for cover content. We would like to thank the reviewers who have taken part in the peerreview process. © 2019 Strojniški vestnik - Journal of Mechanical Engineering. All rights reserved. SV-JME is indexed / abstracted in: SCI-Expanded, Compendex, Inspec, ProQuest-CSA, SCOPUS, TEMA. The list of the remaining bases, in which SV-JME is indexed, is available on the website. The journal is subsidized by Slovenian Research Agency. Strojniški vestnik - Journal of Mechanical Engineering is available on https://www.sv-jme.eu. Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4 Contents Contents Strojniški vestnik - Journal of Mechanical Engineering volume 65, (2019), number 4 Ljubljana, April 2019 ISSN 0039-2480 Published monthly Papers Andrej Jeromen, Edvard Govekar: Nonlinear Dynamic Force Balance Mass-spring-damper Model of Laser Droplet Generation from a Metal Wire 201 Mehmet Direk, Mehmet Selguk Mert, Eren Soylu, Fikret Yüksel: Experimental Investigation of an Automotive Air Conditioning System Using R444A and R152a Refrigerants as Alternatives of R134a 212 Youhang Zhou, Yong Li, Hanjiang Liu: Feature Enhancement Method for Drilling Vibration Signals by Using Wavelet Packet Multi-band Spectral Subtraction 219 Alexandra Catalin Filip, Cristin Olimpiu Morariu, Laurentiu Aurel Mihail, Gheorghe Oancea: Research on the Surface Roughness of Hardox Steel Parts Machined with an Abrasive Waterjet 230 Karolis Januševičius, Juozas Bielskus, Vytautas Martinaitis: Functionality Assessment of Building a Microclimate System Utilising Solar Energy in a Cold Climate 238 Kaikai Luo, Yong Wang, Houlin Liu, Matevž Dular, Jie Chen, Zilong Zhang: Effect of Coating Thickness on a Solid-Liquid Two-Phase Flow Centrifugal Pump under Water Medium 251 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 © 2019 Journal of Mechanical Engineering. All rights reserved. D0l:10.5545/sv-jme.2019.6062 Original Scientific Paper Received for review: 2018-09-03 Received revised form: 2019-02-15 Accepted for publication: 2019-02-21 Nonlinear Dynamic Force Balance Mass-spring-damper Model of Laser Droplet Generation from a Metal Wire Andrej Jeromen* - Edvard Govekar University of Ljubljana, Faculty of Mechanical Engineering, Slovenia In the laser droplet generation process a metal wire is continuously fed into the focus of a pulsed annular laser beam that melts the wire-end and forms a growing pendant droplet that is then detached from the wire. The process is highly complex due to its non-stationarity, non-linearity, and the interplay of numerous physical phenomena. With the aim of describing the process, a low-dimensional, non-linear, dynamic force balance, mass-spring-damper model of the pendant droplet with time-dependent coefficients was formulated based on experimental observations of the process. A comparison between the modelled and experimental droplet centroid vertical position time series and their time-frequency maps showed that the model captures the essential pendant-droplet dynamics in the selected laser-pulse frequency range between 60 Hz and 190 Hz. It was also found that the modelled time of detachment and the detached-droplet diameter were in good agreement with the experimental results, including the bifurcation at the laser-pulse frequency of 120 Hz and the coexistence of two detached-droplet diameter values below that frequency. In addition, the pendent droplet's lateral oscillation and the Rayleigh-Plateau instability were identified as having a significant influence on the process outcome in certain laser-pulse frequency ranges. Keywords: laser droplet generation, mass-spring-damper system, nonlinear model, time series, time-frequency analysis Highlights • Low-dimensional model captures the essential dynamics of the droplet-generation process. • Modelled time series and their time-frequency maps correspond well with experiments. • Model predicts the detached droplet diameter as a function of laser-pulse frequency. • Lateral oscillation and Rayleigh-Plateau instability can influence the process. 0 INTRODUCTION Direct laser metal deposition is one of the established metal-based additive manufacturing processes [1]. It usually employs input material in shape of wire [2] or powder [3]. In contrast to the wire or powder deposition process, where the laser is used to generate melt pool on the workpiece surface, the process of metal droplet deposition significantly reduces the workpiece heat load, since the laser beam is used only for droplet generation. In consequence, metal droplets have been employed in various innovative technologies, including the droplet joining of electrical contacts, and of dissimilar and temperature-sensitive materials [4], as well as 3D printing [5]. A necessary basis for droplet-based technologies is a reliable method for the generation of metal droplets. To achieve this, a number of metal-droplet generation methods based on melt ejection from a nozzle have been developed, employing different techniques to generate the ejection pulse: piezoelectric [6], pneumatic [7], solenoid [8], and electromagnetic [9]. The above methods employ a heated crucible to maintain a reservoir of molten metal. Metals with relatively low melting temperatures are used, such as solder [9], indium, lead, zinc [7], and aluminum [10]. To avoid an energy-inefficient crucible and simplify the system, a method employing induction melting of the wire tip in the pressurized tube leading to the nozzle was proposed [11]. Common to all of the above methods is the nozzle, the diameter of which should be well below the diameter of the generated droplets and is subject to wear. To avoid melt oxidation, leading to nozzle clogging, a low-oxygen-content atmosphere should be provided [7]. Besides the methods based on melt ejection, several laser-based methods have been developed that take advantage of the high level of spatial and temporal control of the laser's energy delivery. They differ in terms of the shape of the input metal material, which can be either a foil [12], a spherical preform [13], or a wire [14]. Although it achieves lower droplet-generation frequencies compared to melt-ejection methods, the laser droplet generation from a wire using an annular laser beam [15] shows strong application potential, since it can employ metals with a high melting point, such as nickel [14], is not very sensitive to oxygen content, and the related system includes no wearable parts. In annular laser beam droplet generation (ALDG), a metal wire is fed axially into the center of the annular laser beam, which is focused on the circumference of the wire-end. ALDG can be performed either as drop-on-demand ALDG [14], where the pendant-droplet formation and detachment *Corr. Author's Address: University of Ljubljana, Faculty of Mechanical Engineering, Aškerčeva 6, SI-1000 Ljubljana, Slovenia, andrej.jeromen@fs.uni-lj.si 201 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 are achieved by a single time-dependent, laser-pulse synchronized with the feeding of the wire, or as continuous ALDG [15], where periodic laser pulses melt the end of the wire, which is fed continuously with a constant velocity into the focus of the annular laser beam. An experimental study of the continuous ALDG has shown that, depending on the laser-pulse frequency, the detachment of the pendant droplet is stimulated by the resonance of a droplet-oscillation eigenmode with the polar wavenumber n = 1, 2, or 3 [15]. In addition, a recent experimental study, based on time series and focused on the n = 1 droplet-oscillation eigenmode in the laser-pulse frequency range between 60 Hz and 190 Hz [16], showed that the essential dynamics of the pendant droplet are very similar to that of a time-varied, mass-spring-damper system. In the study, the effective spring stiffness was determined as a function of the droplet volume and the resonance-caused droplet detachment was confirmed. These results indicate that despite the process being highly complex due to its non-stationarity, non-linearity, and the interplay of numerous physical phenomena, a low-dimensional model could describe the core dynamics of the process. Such a model would provide insight and an understanding of the process, and would also be beneficial for process monitoring, forecasting, and control. Low-dimensional, mass-spring-damper models were already employed for the description of a process similar to continuous ALDG, i.e., gas metal arc welding (GMAW) in the globular mode. The structure of the employed models ranges from the simplest static force balance model that uses either a point mass and constant coefficients [17] (or also includes droplet geometry and time-dependent coefficients [18]) to a dynamic force balance model that includes the inertial force in the droplet-detachment criterion [19]. The most important coefficients of those models are the nonlinear spring stiffness and the damping coefficient, which could be theoretically determined based on the droplet's surface-energy gradient and the viscous damping, respectively. However, the resulting theoretical value of the spring stiffness was found to be too large, which is consistent with the results of the ALDG study [16]. In contrast, the theoretically determined viscous damping coefficient was reported to be much too small to match the experimental droplet dynamics, which could be caused by the rotational fluid flow inside the droplet [18]. Although the above models can provide the basis for formulating the model of continuous ALDG, they cannot be used directly, since the electric arc in GMAW exerts an electromagnetic force that is not present in ALDG as well as there being a difference in the droplet excitation, which in GMAW is achieved by modulation of the electric current, while in continuous ALDG it is achieved by melting of the wire by the laser beam. The aim of this work is to formulate a mathematical model of the continuous ALDG process in the form of a low-dimensional, dynamic force balance, mass-spring-damper system with time-varied coefficients and to evaluate the formulated model by comparing the resulting time series and their time-frequency maps with the experimental ones. To achieve this, the first section briefly describes the continuous ALDG experimental system and process, as well as the relevant results of the previous experimental analysis [16]. Based on these results, in the second section the time-varied coefficients are determined and the dynamic force balance, mass-spring-damper model of droplet growth and detachment is formulated. In the third section, the model-generated time series, their time-frequency maps, and the diameters of the detached droplets are presented, discussed, and compared with the experimental results. 1 EXPERIMENT A schematic presentation of the ALDG experimental system [14] is shown in Fig. 1. The system consists of a Nd:YAG pulsed-laser source, an optical part that transforms the collimated laser beam into an annular form and then focuses it, and a wire-feeding unit, which feeds the wire into the focus of the laser beam. In order to prevent oxidation, the process region is flooded by argon gas of purity 5.0 using a coaxial nozzle. In the continuous ALDG process [15], the wire is fed at a constant feeding velocity, whereas the laser source emits periodic laser pulses at a certain preset frequency in order to melt the wire-end. Due to the periodic heating and melting of the wire directly above the pendant droplet, the droplet grows by incorporating the melted material and is at the same time excited into oscillation since its attachment point is periodically displaced. The continuous ALDG experiments [15] were performed using a nickel wire and the values of the process parameters, summarized in Table 1. At different values of the laser-pulse frequency fp, the power of the laser pulse was adjusted to maintain a constant average laser power Pavg. The process was monitored by means of a high-speed IR sensitive camera, which captured images of the laser beam's focus region with a frequency of 1445 frames per 202 Jeromen, A. - Govekar, E. Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 second, an integration time of 40 ^s, and an image array size of 64x128. Using the captured images, the shape and the position time dependence of the hot pendant droplet were obtained. Fig. 1. Scheme of the ALDG experimental setup Table 1. Summary of the process parameters Parameter Value Wire diameter, dw 0.25 mm Wire feeding velocity, vw 0.06 m/s Laser pulse frequencyf 50 Hz to 300 Hz Laser pulse duration, t 0.8 ms Average laser power, Pavg 120 W resonance of the n = 1 eigenmode [15]. In order to provide the reader with the experimental background that forms the starting point of the present study, the main results of the continuous ALDG experimental study [16] are reproduced and briefly commented on here. In the selected frequency range, the main dynamics of the pendant droplet could be described as that of a forced, time-variable, lightly damped, massspring oscillating system, the mass of which increases with time [16]. The effective spring constant k of the system was found to show a power-law dependence on the volume V of the pendant droplet: k (V )=k / \a V V v o y (1) where k0 = 0.011 N/m, V0 = 1 m3, and a = -0.29, whereas the system's damping ratio could not be estimated from the experimental data [16]. The results of the experimental analysis showed that the pendant droplets detach when in resonance with the laser-pulse frequency fp. However, at frequencies f between 60 Hz and 120 Hz smaller droplets were also observed to detach when in resonance with the laser-pulse frequency's second multiple 2fp. In Fig. 2, the detached droplet diameter dd values [15] are reproduced in the selected laser-pulse frequency fp range. The upper branch in Fig. 2 denotes the diameters of the droplets that were detached in resonance with fp, whereas the lower branch denotes the diameters of the smaller droplets that were detached in resonance with 2fp. In the frequency range f between 60 Hz and 120 Hz droplets belonging to both branches were generated in a seemingly random sequence. Fig. 2. Experimental detached-droplet diameter dd as a function of the laser-pulse frequencyfp The present study is focused on the laser-pulse frequency f range between 60 Hz and 190 Hz, where the pendant droplet's detachment is triggered by the 2 MODEL FORMULATION Before formulating the model, several relevant properties of the continuous ALDG process are presented. The pendant droplet, which is going to be represented by a point-mass in the model, grows and oscillates in the process. Due to its increasing mass, the properties of the droplet, such as its eigenmode frequencies, the attachment-to-droplet-diameter ratio, the effective spring constant, and the damping coefficient, are all time dependent. Additionally, due to the periodic heating and melting of the wire above the droplet, the pendant droplet's oscillation is forced by the periodic mass increases and the corresponding attachment-point displacements. It should be noted that the forcing in continuous ALDG is not harmonic, but closer to impulse-like, since the duration of the laser pulse t is short compared to the pulse period 1f between 5.3 ms and 17 ms in the selected frequency Nonlinear Dynamic Force Balance Mass-spring-damper Model of Laser Droplet Generation from a Metal Wire 203 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 range. The pendant droplet can also be described by the values of the following governing dimensionless groups [20], calculated using the properties of molten nickel, summarized in Table 2, and a typical pendant-droplet diameter d = 1 mm. The large value of the Reynolds number Re = (adp/(2j2))1/2 = 1.6x104 indicates that the influence of viscosity is small, while the small value of the gravitational Bond number is G = (pgd2)/(4a) = 1.1x10-2, where g = 9.81 m/s2 is the acceleration due to gravity, indicates that the influence of surface tension is large compared to gravity, which is reflected in the observed, nearly spherical, shape of the pendant droplet. Table 2. Summary of physical properties of molten nickel at the melting temperature [21] The equation of motion of the point mass position y was formulated as a second-order differential equation with time-varying coefficients, as follows: Property Value Temperature of melting, Tm 1728 K Density, p 7905 kg/m3 Viscosity, j 1.66x10-4 Pas Surface tension, a 1.778 N/m Thermal diffusivity, a 5.8X10-6 m2/s Fig. 3. Scheme of the mass-spring-damper model and its geometrical relation to the pendant droplet In order to model the pendant droplet's centroid oscillation and detachment in continuous ALDG, a forced, time-variable, mass-spring-damper model was formulated. The scheme of the model system and its geometrical relation to the pendant droplet is presented in Fig. 3. The model system consists of a time-dependent point mass m(t), which is suspended on a spring with a time-dependent spring stiffness k(t) and an unloaded length l0(t), and a viscous damper with a time-dependent damping coefficient c(t). The positiony(t) of the model point mass m corresponds to the position of the pendant-droplet centroid, whereas the position yj(t) of the model spring-attachment point corresponds to the position of the pendant droplet neck, i.e., of the attachment to the wire. m (t ) y = -m (t ) g + k (t )( (t )-y -10 (t ))--c(t)(yi -y). (2) The terms on the right-hand side of Eq. (2) represent the gravitational force, the spring force, and the viscous damping force. The drag force on the droplet moving in the surrounding protective atmosphere of argon can be taken as negligible. For this reason the damping force was assumed to be mainly the result of the droplet's inner flows [18], and was modelled as being proportional to the relative velocity y1 - y between the droplet centroid and the droplet neck, i.e., the spring attachment position. In the following text, the time-dependent coefficients of Eq. (2), the initial values of the model, and the pendant-droplet detachment condition are determined based on the experimental observations. The displacement of the spring-attachment position y1 is the source of the external forcing. Its time dependence yj(/) is a consequence of the constant-velocity wire feeding combined with the laser-induced periodic melting of the wire above the droplet neck. The function yj(/) was formulated as a piecewise linear function. Between the consecutive pulses, the spring-attachment position yjwas assumed to decrease linearly with the velocity of the wire feeding vw. At the time of the laser-pulse trigger, the laser beam is directed onto the wire circumference above the pendant droplet's neck. For this reason, during a laser pulse of duration t, the position yj was assumed to increase linearly up to the position of the laser beam focus at y=0, as indicated in Fig. 3. The minimum value yjmin of the spring-attachment position yj was thus expressed as yjmin=-vw(t0- t), where vw is the downward wire-feeding velocity magnitude, t0 = 1/fp denotes the time period of the laser pulses, and t is the duration of the laser pulse. The experimentally determined minimum values of y1 correspond well to the above-determined y1min for very small droplets, i.e., at the very beginning of their growth. However, during the time t1 = 0.3 s from the beginning of the droplet growth, the experimentally observed minimum value gradually increased to approximately half of its theoretical value -vw(t0 - t). The observed increase of y1min can be explained by the increasing mass of the liquid pendant droplet, which contains an increased amount of energy and is increasingly able to melt the wire between the laser pulses. However, since the temperature of the droplet is close to the 204 Jeromen, A. - Govekar, E. Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 melting temperature, the melting of the wire by the pendant droplet saturates and the minimum value y1min stabilizes at half of its theoretical value for the selected values of the average laser power, the wire diameter dw, and the feeding velocity vw. In order to take into account this experimentally observed effect, the minimum value y1min was defined as time dependent: >1min (t ) = -Vw (t0 -T) 1-- 2t ; t < ti, .((0-t). 2 ; t > ty (3) Based on this, the spring-attachment position's time dependenceyx(t) was then expressed as: *(t ) = t — nt y 1min (t )"-^ ^ nt0 - t <(" + 1)t0 >1mm (t) 1 - t -(nto-T) ; nt0-t< t <(n + l)t0, (4) where n = 0, 1,2, ... corresponds to the consecutive laser-pulse period. In Fig. 4 the time dependencies of the spring-attachment position y1 and its minimum value y1min are represented by solid and dashed curves, respectively. To avoid non-realistic discontinuities of the spring-attachment position's time derivative y1, the above-defined y1(t) time dependence was smoothed in order to take into account the rise time of the laser source rls = 0.2 ms and the thermal penetration time Ttp, estimated by rtp = (dw/2)2 / (na)) = 0.85 ms. For smoothing, a moving average filter was employed with a length of 0.135 t0, which is of the order of the thermal penetration time Ttp and resulted in the best agreement between the results of the model and the experiment. Fig. 4. Schematic presentation of the spring-attachment position's time dependence y1(t), the full curve; and its minimum value time dependence y1min(t), the dashed curve The pendant droplet's mass time dependence m(t) follows from the melting of the wire, as the melted material is incorporated into the pendant droplet. In the model, the wire is considered to melt whenever the spring-attachment position y1, i.e., the droplet neck, travels downwards more slowly than the wire, which is fed by the constant feeding velocity vw. In this way, in addition to the melting of the wire by the laser pulse, melting of the wire between the laser pulses due to the droplet's internal energy is also taken into account. Thus, during an infinitesimally small time interval dt, the melted wire length dlw can be expressed by a time derivative of the spring-attachment position y1 : dlw = (yy1 + vw)dt. The corresponding pendant droplet's mass increase is then: dm = Aw p dlw = ^wp( y1 + vw)dt, where = ndw2 / 4 is the cross-sectional area of the wire, and p is the density of the wire. The expression for dm is integrated in order to determine the pendant droplet's mass time dependence m(t): t m () = m0 + j^w + Vw (5) 0 where m0 = m(t = 0) is the pendant droplet's initial mass. Fig. 5 provides a schematic representation of the pendant droplet's mass time dependence m(t), which corresponds to the spring-attachment position's time dependence y1(t), as presented in Fig. 4. From Fig. 5 it can be seen that at the beginning the mass of the pendant droplet increases only during the laser pulse of duration t. However, at later times, the mass slowly increases between the laser pulses due to melting of the wire by a large and hot pendant droplet. At all times, the droplet mass increase Am during one laserpulse period is equal to the mass of the wire that is fed during the same time: Am = Aw p vw t0, as indicated in Fig. 5. Fig. 5. Schematic presentation of pendant droplet's mass time dependence m(t), corresponding to the spring-attachment position's time dependencey1 (t) The length l0 of the unloaded spring was presumed to be equal to the radius d/2 of a spherical t Nonlinear Dynamic Force Balance Mass-spring-damper Model of Laser Droplet Generation from a Metal Wire 205 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 droplet, which is in accordance with the small gravitational Bond number G and the experimental observations. The time dependence of the length l0 was determined using the time dependence of the radius d(t)/2, expressed by the droplet mass m(t) and its density p: d (t lo (t) = ■ 2 3m (t ) (6) The time-dependent spring stiffness k(t) was determined from the experimental effective spring stiffness volume dependence k(V) in Eq. (1) [16], by expressing the droplet volume V(t) by its mass m(t) and density p: k (t ) = k0 m (t ) ~PV~o (7) Since in previous work the value of the viscous damping coefficient was found to be too small to describe the damping of the metal pendant droplet, presumably due to rotational flow in the droplet, and, therefore, a constant value damping coefficient was used [18], the damping coefficient c(t) of the proposed model was assumed to be constant. The value of 1.5*10-3 kg/s used in the calculation most closely reproduced the experimental droplet-oscillation amplitudes and detached-droplet diameters, and is similar to the value of 2.0*10-3 kg/s, which was used in a similar model of a pendant steel droplet [18]. The employed value is indeed much larger than the viscous damping coefficient, estimated using the expression / (d/2) [19], which for a nickel droplet of diameter d = 1 mm gives 5.0*10-5 kg/s. This indicates the important influence of other possible factors, such as rotational flow in the oscillating pendant droplet, driven by the temperature gradient, which are beyond the scope of this study. The initial pendant droplet's mass m(t = 0) = m0 corresponds to the mass of the remaining pendant part after the detachment of the preceding droplet. The mass of the remaining pendant part was experimentally hard to estimate from the recorded infrared (IR) intensity images due to its small diameter and low temperature. For this reason, an initial pendant-droplet mass m0 value was selected in the interval of less than 2.4*10-7 kg, which was experimentally estimated as the upper initial mass limit, and is of the order of the typical droplet massincrease Am during one laserpulse period. In order to define the initial state of the model, the following initial values were selected. The initial spring-attachment position yx(t = 0) was chosen to be 0. The initial point-mass position y(t = 0) was then determined as its static force balance position: y(t = 0) = l0(t = 0) - m0 g / k(t=0). The initial point-mass velocity y was defined as being equal to the negative value of the wire-feeding velocity magnitude y (t = 0) = -vw. In this way, at t=0 the pendant droplet moves downwards with the wire. The detachment of the pendant droplet was modelled using a dynamic force balance condition. For pendant droplet's detachment, the ratio n of the combined inertial, my, and gravitational, mg, forces to the maximum surface-tension force, ndwa, which retains the pendant droplet attached to the wire, should exceed a value of one. Additionally, the point-mass velocity y should at the same time be directed downwards in order to avoid the case of detachment and immediate re-attachment of an upwards travelling droplet, which was observed experimentally [15]. Thus, the modelled pendant droplet detaches when both of the following conditions are fulfilled: > 1, nd a y < 0. 3 RESULTS AND DISCUSSION (8) The formulated model was solved using the Runge-Kutta numerical integration method at laser-pulse frequency f values between 60 Hz and 190 Hz, using a step of 10 Hz. At each frequency value, the model was solved for 31 equidistant values of the initial pendant-droplet mass m0 from 0.04*10-7 kg up to 2.4*10-7 kg, where the latter value was experimentally estimated as the upper initial mass limit. In Fig. 6, examples of the modelled point-mass vertical position time series y(t) are presented as thick red lines for laser-pulse frequency fp values of (70, 100, 140, and 180) Hz. The initial pendant-droplet mass m0 was, in all cases, 0.748*10-7 kg. For comparison, the examples of the corresponding experimental pendant-droplet centroid vertical position time series y(t) estimated from the high-speed IR records of the continuous ALDG process [16] are plotted in Fig. 6 using thin blue lines. In the presented cases the droplet detached is in resonance with fp. The moment of the detachment is denoted by a red square and a blue circle for the model and experiment, respectively. As can be seen in Fig. 6, the model successfully reproduces the experimental vertical centroid position time series y(t) except in some sections of a 206 Jeromen, A. - Govekar, E. Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 the droplet growth, where an obvious discrepancy can be observed. The cause of this discrepancy is the pendulum-like lateral oscillation of the pendant droplet, which is revealed by the experimental horizontal centroid position time series x(t), denoted by a dotted black line in Fig. 6. Indeed, at the time of the observed discrepancy, a large amplitude of lateral oscillation can be seen in the horizontal centroid position time series x(t). Here, it should be noted that only a projection of the lateral oscillation on the IR camera image plane could be observed. In the case of a lateral oscillation perpendicular to the image plane, an oscillation in the x coordinate would not be detected, while the lateral oscillation would still influence the detected vertical oscillation. The observed lateral oscillation could be induced by a small deviation from the axial symmetry, e.g., a small misalignment between the wire and the laser beam, and could result in large amplitudes in the resonance of one of the many oscillation modes of the pendant-droplet system, which is very similar to a spring pendulum [22]. Since the presented one-dimensional, mass-spring-damper model includes only the vertical oscillation, the lateral degree of freedom is unaccounted for, which probably led to the above-described discrepancy. The modelled detachment time falls within 5 % of the experimental detachment time, except at laser-pulse frequencies f above 150 Hz, where the modelled detachment time shown in Fig. 6, f = 180 Hz, is much too small, which will be commented on at the end of the section. The good agreement between the modelled and experimental times of detachment is a consequence of the pendant-droplet resonance with the laser-pulse frequency fp, which provides both a large enough oscillation amplitude and a favorable phase difference of n/2 between the pendant-droplet oscillation and the laser pulses [16]. The phenomenon is illustrated in Fig. 7, Fig. 6. Examples of the vertical centroid position time series y(t) resulting from the model (marked by thick red lines, with small squares denoting the moment of detachment), and from the experiment (marked by thin blue lines, with circles denoting the moment of detachment) with the experimentally determined horizontal centroid position time series x(t) (marked by dotted black lines) of the droplets that detached in resonance with fp Nonlinear Dynamic Force Balance Mass-spring-damper Model of Laser Droplet Generation from a Metal Wire 207 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 where the harmonic-wavelet-based, time-frequency map [16] of the vertical position time series y(t) is presented for the laser-pulse frequencies fp = 70 Hz and 180 Hz for both the modelled case (Fig. 7a)and the experimental case (Fig. 7b). In this figure, the small amplitudes are denoted by dark tones, whereas the large amplitudes are denoted by bright tones. The modelled and experimental cases presented in Fig. 7 show the same time evolution of the resonances at lower and lower multiples of the laser-pulse frequency fp, which is ended by detachment, denoted by a vertical dotted line, triggered by the resonance at the laser-pulse frequency fp. The detachment in all the presented cases takes place after the onset of the peak of the vertical oscillation amplitude, because only then does the phase difference between the droplet oscillation and the laser pulses reach a value close to nil [16], when the inertial force is large enough and the detachment conditions are fulfilled. The non-monotonically increasing amplitude at f=fp, which can be observed in the experimental time-frequency maps in Fig. 7b, follows from the vertical position y(t) amplitude variation due to the lateral oscillation. In the experiment as well as in the model, some smaller droplets were observed to detach at laser-pulse frequency fp values below 1l0 Hz. The detachment of these droplets was shown to coincide with the vertical oscillation resonance corresponding to the second multiple of the laser-pulse frequency f [16]. For droplets detached in this way at laser-pulse frequency fp values of 70 Hz and 100 Hz, examples of the model point-mass vertical position time series y(t) are presented in Fig. 8 by means of a thick red line. The initial pendant-droplet mass m0 was, in both cases, 0.512*10-7 kg. For comparison, the corresponding experimental pendant-droplet centroid vertical y(t) and the horizontal x(t) position time series are plotted by means of a thin blue line and a dotted black line, respectively. The moment of detachment is denoted by a small red square for the model and a blue circle for the experiment. Again, the model successfully reproduced the experimental vertical centroid position time series y(t), where the horizontal oscillation amplitude is low. It can be seen from Fig. 8 that, before detachment, as well as at the moment of detachment, the horizontal oscillation amplitude is high. It has been experimentally observed that the triggering of the detachment near the lfp resonance is due to the resonance of the oscillation mode combining vertical and lateral degrees of freedom [16]. Thus, the modelled time of detachment, where only the vertical oscillation is taken into account, is less similar to the experimental time of detachment. Also, the comparison of the corresponding time-frequency plots for the modelled case (Fig. 9a) and the experimental Fig. 7. Time-frequency maps of a) the modelled, and b) the experimental vertical centroid position time series y(t) atfp-=70 Hz and 180 Hz of the droplets that detached in resonance with fp 208 Jeromen, A. - Govekar, E. Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 Fig. 8. Examples of vertical centroid position time seriesy(t) resulting from the model (marked by a thick red line, with small squares denoting detachment) and the experiment (marked by a thin blue line, with small circles denoting the detachment) with the experimental horizontal centroid position time series x(t) (shown by a dotted black line) of the droplets that detached in resonance with 2fp Fig. 9. Time-frequency maps of a) the modelled and b) the experimental vertical centroid position time series y(t) at fp =70 Hz of the droplets that detached in resonance with 2f case (Fig. 9b) at the laser-pulse frequency fp = 70 Hz reveals this difference. Whereas in the modelled case (Fig. 9a) the detachment, denoted by a vertical dotted line, takes place immediately after the resonant peak at 2fp = 140 Hz, when the phase difference between the droplet oscillation and the laser pulses reaches a value close to n/2, in the experimental case (Fig. 9b), the same peak is less pronounced, and the detachment is triggered later since it is a consequence of a resonance of the combined oscillation mode that includes both vertical and lateral degrees of freedom. From the point of view of applications, the diameter dd of the detached droplet is one of the most important output parameters of the droplet-generation process. Since the pendant-droplet detachment in continuous ALDG at laser-pulse frequency fp values between 60 Hz and 190 Hz is triggered by resonance, the detached droplet diameter dd depends on the laser-pulse frequency fp. In order to show this dependence, the modelled and experimental detached-droplet diameters dd are presented in Fig. 10 as small red circles and blue points, respectively, with their dependence on the laser-pulse frequency fp. As is evident from the figure, in the laser-pulse frequency fp interval between 60 Hz and 150 Hz, the model successfully reproduces the experimentally observed diameters, corresponding to the detachment triggered by the resonance at the laser-pulse frequency fp, denoted by the higher branch of blue points, and by the resonance at the second multiple 2fp, denoted by the lower branch of blue points. Also, the bifurcation near 120 Hz, below which the two branches coexist, was reproduced by the model. In Fig. 10, the largest discrepancy between the modelled and the experimental values of the detached-droplet diameter dd can be observed at the lowest and at high values of the laser-pulse frequency fp. At the lowest frequency fp = 60 Hz, the model predicts Nonlinear Dynamic Force Balance Mass-spring-damper Model of Laser Droplet Generation from a Metal Wire 209 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 droplets with a diameter dd between 1.2 and 1.3 mm, whereas in the experiment the droplet diameter dd shows much larger variations. The reason for this discrepancy is the influence of the Rayleigh-Plateau instability of the molten column of wire above the droplet, which at this low value of the laserpulse frequency fp is long enough directly after the laser pulse to become unstable and can cause the detachment of a droplet with a smaller diameter [15]. On the other hand, at high values of the laser-pulse frequency, fp > 150 Hz, the modelled values of the detached droplet diameter dd are clearly smaller than the experimentally observed values, which is also evident in Fig. 6 in the case of fp = 180 Hz. Due to the high laser-pulse frequency, the vertical forcing peak-to-peak amplitude LVimini is smaller, and the pendant-droplet detachment triggered by a vertical oscillation can be delayed due to the onset of the lateral oscillation, which slows down the rate of increase of the vertical oscillation amplitude. As can be seen in Fig. 10, up to the laser-pulse frequency ./¿=180 Hz, some experimental droplets were observed to detach with approximately the same diameter dd (and at the same time) as the modelled ones. In these cases a lateral oscillation was not present. At the modelled laser-pulse frequency values f > 170 Hz, the model predicts detached droplets with a diameter dd near 1.0 mm, whereas, in the experiment, also droplets with diameters greater than 2.0 mm were observed [15]. The reason for this discrepancy is the low vertical forcing amplitude combined with the lateral oscillation of the pendant droplet, which enable the pendant droplet to pass the resonance at fp without detachment, so that it detaches later, having a larger diameter, due to the resonance of the n = 2 mode [15]. Fig. 10. Diameter dd of the modelled and experimentally observed detached droplet In dependence on the laser pulse frequency fp 4 CONCLUSIONS In the paper, the pendent-droplet oscillation and detachment in the continuous ALDG process were described by means of a mathematical, dynamic force balance, mass-spring-damper model in the laser-pulse frequency range between 60 Hz and 190 Hz. The time-dependent coefficients of the model were determined based on the experimental observations. The comparison between the modelled and experimental droplet centroid vertical position time series and their time-frequency maps show that the proposed model was able to reproduce the essential dynamics. The model also reproduced the pendant-droplet growth and detachment, triggered by the n = 1 eigenmode resonance of the pendant droplet with the frequency of the laser pulses or its second multiple. Compared to the experimental results, the detached droplet diameter as a function of the laser-pulse frequency was predicted well by the model, including the bifurcation at 120 Hz and the coexistence of two detached-droplet diameter values at laser-pulse frequencies below 120 Hz. However, the modeled detached-droplet diameter differs from the experimental values at both ends of the selected laser-pulse frequency interval, where the droplet detachment is influenced by phenomena that are not accounted for by the model, but can be explained based on the experimental observations. This is because, at the lowest frequency, the droplet detachment is influenced by the Rayleigh-Plateau instability, whereas at the highest frequencies, the onset of the lateral pendulum-like oscillation can delay or even prevent the n = 1 eigenmode resonant detachment. The description of the continuous ALDG could be improved by including the lateral degree of freedom in the presented model. 5 ACKNOWLEDGEMENTS The authors acknowledge the financial support provided by the Slovenian Research Agency (Research core funding No. P2-0241). 6 REFERENCES [1] Schmidt, M., Merklein, M., Bourell, D., Dimitrov, D., Hausotte, T., Wegener, K., Overmeyer, L., Vollertsen, F., Levy, G. N. (2017). Laser based additive manufacturing in industry and academia. CIRP Annals, vol. 66, no. 2, p. 561-583, D0l:10.1016/j.cirp.2017.05.011. [2] Kotar, M., Govekar, E. (2018). The influence of the workpiece illumination proportion in annular laser beam wire deposition process. Procedla CIRP, vol. 74, p. 228-232, D0I:10.1016/j. procir.2018.08.100. 210 Jeromen, A. - Govekar, E. Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 201-211 [3] Mutiu, E.F., Akinlabi, E.T. (2018). Influence of laser power on improving the wear properties of laser-deposited Ti-6Al-4V+B4C composite. Strojniški vestnik - Journal of Mechanical Engineering, vol. 64, no. 7-8, p. 488-495, D0I:10.5545/sv-jme.2018.5362. [4] Govekar, E., Jerič, A., Weigl, M., Schmidt, M. (2009). Laser droplet generation: Application to droplet joining. CIRP Annals - Manufacturing Technology, vol. 58, no. 1, p. 205-208, D0I:10.1016/j.cirp.2009.03.005. [5] Chen, C.-A., Chun, J.-H., Sohlenius, G. (1997). Development of a droplet-based manufacturing process for free-form fabrication. CIRP Annals, vol. 46, no. 1, p. 131-134, D0I:10.1016/S0007-8506(07)60791-4. [6] Lee, T.-M., Kang, T.G., Yang, J.-S., Jo, J., Kim, K.-Y., Choi, B.-O., Kim, D.-S. (2008). Drop-on-demand solder droplet jetting system for fabricating microstructure. IEEE Transactions on Electronics Packaging Manufacturing, vol. 31, no. 3, p. 202210, D0I:10.1109/TEPM.2008.926285. [7] Cheng, S.X., Li, T., Chandra, S. (2005). Producing molten metal droplets with a pneumatic droplet-on-demand generator. Journal of Materials Processing Technology, vol. 159, no. 3, p. 295-302, D0I:10.1016/j.jmatprotec.2004.05.016. [8] Sohn, H., Yang, D.Y. (2005). Drop-on-demand deposition of superheated metal droplets for selective infiltration manufacturing. Materials Science and Engineering: A, vol. 392, no. 1-2, p. 415-421, D0I:10.1016/j.msea.2004.09.049. [9] Luo, Z., Wang, X., Wang, L., Sun, D., Li, Z. (2017). Drop-on-demand electromagnetic printing of metallic droplets. Materials Letters, vol. 188, p. 184-187, D0I:10.1016/j. matlet.2016.11.021. [10] Yi, H., Qi, L., Luo, J., Zhang, D., Li, H., Hou, X. (2018). Effect of the surface morphology of solidified droplet on remelting between neighboring aluminum droplets. International Journal of Machine Tools and Manufacture, vol. 130-131, p. 1-11, D0I:10.1016/j.ijmachtools.2018.03.006. [11] Vega, E.J., Cabezas, M.G., Muñoz-Sánchez, B.N., Montanero, J.M., Gañán-Calvo, A.M. (2014). A novel technique to produce metallic microdrops for additive manufacturing. The International Journal of Advanced Manufacturing Technology, vol. 70, no. 5-8, p. 1395-1402, D0I:10.1007/s00170-013-5357-3. [12] Jeromen, A., Kuznetsov, A., Govekar, E. (2014). Laser droplet generation from a metal foil. Physics Procedia, vol. 56, p. 720729, D0I:10.1016/j.phpro.2014.08.079. [13] Stein, S., Dobler, M., Radel, T., Strauß, M., Breitschwerdt, H., Hugger, F., Roth, S., Schmidt, M. (2017). Experimental and numerical investigations regarding laser drop on demand jetting of Cu alloys. Production Engineering, vol. 11, no. 3, p. 275-284, D0I:10.1007/s11740-017-0737-4. [14] Govekar, E., Kuznetsov, A., Jerič, A. (2016). Drop on demand generation from a metal wire by means of an annular laser beam. Journal of Materials Processing Technology, vol. 227, p. 59-70, D0I:10.1016/j.jmatprotec.2015.07.026. [15] Kuznetsov, A., Jeromen, A., Govekar, E. (2014). Droplet detachment regimes in annular laser beam droplet generation from a metal wire. CIRP Annals - Manufacturing Technology, vol. 63, no. 1, p. 225-228, D0I:10.1016/j.cirp.2014.03.051. [16] Jeromen, A., Govekar, E. (2016). Time series analysis based study of a mass-spring-like oscillation and detachment of a metal pendant droplet. Mechanical Systems and Signal Processing, vol. 80, p. 503-516, D0I:10.1016/j. ymssp.2016.04.038. [17] Watkins, A.D., Smartt, H.B., Johnson, J.A. (1992). A dynamic model of droplet growth and detachment in GMAW, Proceedings of 3rd International Conference on Trends in Welding Research, Gatlingburg, p. 993-997. [18] Jones, L.A., Eagar, T.W., Lang, J.H. (1998). A dynamic model of drops detaching from a gas metal arc welding electrode. Journal of Physics D: Applied Physics, vol. 31, no. 1, p. 107123, D0I:10.1088/0022-3727/31/1/014. [19] Choi, J.H., Lee, J., Yoo, C.D. (2001). Dynamic force balance model for metal transfer analysis in arc welding. Journal of Physics D: Applied Physics, vol. 34, no. 17, p. 2658-2664, D0I:10.1088/0022-3727/34/17/313. [20] Basaran, O.A., DePaoli, D.W. (1994). Nonlinear oscillations of pendant drops. Physics of Fluids, vol. 6, no. 9, p. 2923-2943, D0I:10.1063/1.868120. [21] Liu, H. (2000). Science and Engineering of Droplets: Fundamentals and Applications, Noyes Publications, Park Ridge & William Andrew Publishing, Norwich. [22] van der Weele, J.P., de Kleine, E. (1996). The order-chaosorder sequence in the spring pendulum. Physica A: Statistical Mechanics and its Applications, vol. 228, no. 1-4, p. 245-272, D0I:10.1016/0378-4371(95)00426-2. Nonlinear Dynamic Force Balance Mass-spring-damper Model of Laser Droplet Generation from a Metal Wire 211 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 212-218 © 2019 Journal of Mechanical Engineering. All rights reserved. D0l:10.5545/sv-jme.2019.6040 Original Scientific Paper Received for review: 2018-09-03 Received revised form: 2019-02-15 Accepted for publication: 2019-02-21 Experimental Investigation of an Automotive Air Conditioning System Using R444A and R152a Refrigerants as Alternatives of R134a Mehmet Direk* - Mehmet Selguk Mert - Eren Soylu - Fikret Yuksel Yalova University, Faculty of Engineering, Turkey The conventional refrigerants used in refrigeration systems cause global warming effect if released into the atmosphere. To overcome this problem, refrigerants are classified according to their global warming potential (GWP) values, and some restrictions have been placed on the use and transport of these refrigerants with European regulations. Vapor compression refrigeration (VCR) systems are used to provide desired conditions in many systems and as well as in automobiles. Most of the automotive air conditioning (AAC) systems use R134a which has GWP value of 1300. In this study, the use of the refrigerants which have low GWP values, namely, R152a and R444A as an alternative to the R134a for AAC systems, were investigated experimentally. Moreover, the effect of a plate-type internal heat exchanger (IHX) added into the system was studied with the refrigerant R444A. The cooling capacity, compressor power, coefficient of performance (COP) and total exergy destruction have been determined as performance parameters. Based on the obtained results, the highest COP and cooling capacity were achieved when the R152a used as the refrigerant. Furthermore, the COP of the system was increased when the IHX was used, with R444A as the refrigerant. Keywords: R134a, R152a and R444A, automotive air conditioning, coefficient of performance Highlights • The operation of an AAC system using R134a, R152a and R444A was tested. • The influence of IHX on the system performance was revealed with R444A. • The highest COP was obtained when the system was working with R152a. • IHX based system with R444A (with IHX) had the lowest exergetic COP. 0 INTRODUCTION The European Union brings some restrictions to the use of some refrigerants, which has been higher than 150 global warming potential (GWP) value, in mobile air conditioning (MAC) systems as in directive no 2006/40/EC [1]. Then these GWP limits were extended air conditioning and refrigeration systems with regulation the Regulation (EU) No, 517/2014 [2]. Presently, most of the MAC systems use R134a as the refrigerant; however, the GWP value of the R134a is 1300 [3] and [4]. The refrigerants such as R1234yf, R1234ze(E), R152a and R444A, can be considered as alternative refrigerants having low GWP for automotive air conditioning systems. Table 1 shows the thermophysical properties of the above-mentioned refrigerants. Among them, R1234ze(E) can be used instead of R134a due to the their similar properties [5]. However, the cooling performance of R1234ze(E) is 30 % lower in average than that of R134a when used in a similar medium capacity vapor compression system [6] to [8]. Furthermore, the performance of R1234ze(E) can be improved when the mixture of some refrigerants used to obtain required properties. The R444A consists of 83 % R1234ze(E), 12 % R32 and 5 % R152a, (by mass) and its GWP value is 93 [4] and [9] which meets the European regulations. Devecioglu and Orug [10] calculated the performance parameters of R1234yf, R444A and R445A for a MAC system. It was found that the R444A and R445A have lower cooling capacity, but higher coefficient of performance (COP) than that of R1234yf. Lee et al. [11] investigated the performance of R444A, R445A, R152a, and R1234yf in an automotive air conditioning (AAC) system as an alternative to R134a. The highest COP was obtained when the system was working with R152a. Cheng et al. [9] tested an air source heat pump (HP) system using the different concentrations of R32/R1234ze(E) mixture numerically. It was determined that the heating and cooling capacities were improved compared to the R134a baseline system. Meng et al. [12] determined the performance of a refrigeration system using the mixture of R152a and R1234ze(E) refrigerants as an alternative to R134a. The mixture of 50 % R1234ze(E) and 50 % R152a was found to be the best alternative for R134a. It was demonstrated that the cooling capacity of R1234ze(E)/R152a mixture was very similar to R134a and can be used without any change in the compressor. Li et al. [13] reported that R134a 212 *Corr. Author's Address: Yalova University, Faculty of Engineering, Energy Systems Engineering Department, Yalova, Turkey, mehmetdirek@hotmail.com Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 212-218 Table 1. Refrigerant properties of R134a, R152a, R444A, R1234ze(E), R1234yf [3] Property R134a R152a R444A R1234ze(E) R1234yf ASHRAE safety classification_A_A2_A2L_A2L_A2L Ozone depletion potential (ODP) 0 0 0 0 0 GWP [4]_1300_124_93_4__4 Critical temperature [°K]_374.21_386.41_374.39_382.51_367.85 Critical pressure [kPa]_4059.3_4516.8_4235.8_3634.9_3382.2 Vapor density [kg-m-3] (at 0 °C)_32.35_1847_285_26.32_37.92 Liquid density [kg-m-3] (at 0 °C)_1294.8_959.11_1199.1_1240.1_1176.3 Latent heat of vaporization [kJ-kg-1] (at 100 kPa) 177.78 279.36 180.5 166.92 145.37" and R152a work at similar pressure ranges because of the close boiling points (R134a -26.09 °C and R152a -24.05 °C). R152a has 55 % higher latent heat value and 40 % lower liquid density when compared with R134a. Ghodbane [14] tested R152a in an AAC system as an alternative to the R134a and achieved a 15 % higher COP. In another study, Scherer et al. [15] tested an HP system in which R134a and R152a were used as refrigerants and obtained similar results. Lee [16] tested an HP system which was using R134a for zero emission vehicles. The results demonstrated that when the temperature was -10 °C, the heating capacity was insufficient for cabin heating. The performance of a vapor compression refrigeration (VCR) system can be increased by installing an internal heat exchanger (IHX) into the liquid and gas suction line. The IHX increases subcooling degree so that the enthalpy difference in the evaporator is increased and the cooling capacity is increased as well [17]. Moles et al. [18] evaluated the performance of a refrigeration system using R1234ze(E) and R1234yf as an alternative of R134a and found that the cooling capacity increases when the IHX is activated, and the COP was found to be very close to that of R134a. Mota-Babiloni et al. [19] observed that the mass flow rate was 18 % to 29 % lower than that of R134a when the system was working with R1234ze(E). They also determined that R1234ze(E) can provide the same cooling capacity as R134a at 43 % higher compressor speeds with IHX. In another study, a VCR system using IHX with R450A as an alternative to R134a was investigated. It was found that the IHX influence positively to the COP of R450A, in a closed value of R134a [20]. In this study, the low GWP refrigerants, namely, the R152a and R444A, were used as an alternative to R134a, to investigate the performance of an AAC system. An IHX was also used to examine the performance of the system, and to provide an improvement. The experimental AAC system was employed with on and off valves to activate the liquid and gas suction line of the IHX. In the first part, the R134a, R152a and R444A were tested on an AAC system designed for R134a without IHX. Additionally, the effect of the IHX on the system performance was evaluated when the R444A used as the refrigerant. 1 THERMODYNAMIC ANALYSIS The general energy balance for the components in the test system can be written as in Eq. (1). QCV - WCV = £ mh - £ moho . (1) The compressor power and cooling capacity can be evaluated by Eqs. (2) and (3), respectively. Wcomp = mr (h - h), (2) Qevap = mr (( - K). (3) The COP can be determined, as follow: COP = Q / W . (4) z-'evap comp In thermodynamic systems, it is very critical that how much of the energy is converted. When these systems are evaluated in regarding to environmental conditions, this evaluation is expressed by the concept of exergy. The reference states were taken as T0 = 298.15 K and P0 = 1 atm for this study. General exergy destruction equation for closed control volume can be written as follows. ¿xd-Tjq-ztw(5) The exergy destruction rate for each component namely; compressor, thermostatic expansion valve, evaporator, internal heat exchanger was evaluated, as follows: Experimental Investigation of an Automotive Air Conditioning System Using R444A and R152a Refrigerants as Alternatives of R134a 213 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 212-218 Exd ,comp = mr (l - ¥2 ) + Komp , (6) Exd ,TXV = ™r (3 -Va)> (7) Exd ,evap = ™ r (V l -V 4 ) + ™ al (5 - V6 ) , (8) Exd,cond = ™r (3 -¥2 ) + ™a2 (7 -VS), (9) Exd,IHX = ™r (lb ~Vla ) + ™r (3b 3a ). (10) Total exergy destruction can be found from the Eq. (11). Exd,t = Exd,comp + Exd,cond + Exd,TXV ' +Exd evap + Exd iiHX . (11) The exergetic efficiency can be calculated from Eq. (12) [21]. Qe\ i T\ 1 - T T v 1l y -Qco T 1 - T T v h y m w + wv airr air,inlet comp (12) Exergetic COP can be determined from the following equations: T T (13) COP t = carnot th - Tl COP COP __actual Exergetic COP - 2 EXPERIMENTAL SETUP (14) The experimental setup of the AAC system can be seen in Fig. 1. The system also includes a platetype heat exchanger between the liquid and gas suction line. The properties of experimental system components were given in Table 2. iAIRl ® _TLTL © © RESISTANCE HEATER (^FLOWMETER RECEIVER JrV I-1 TANK IHX r • '; '©e * t . ffid -ta© TXV © -Dti— Si CONDENSER <Ê>© » © LH ¿AIR J 1 5 © JTTL © © EVAPORATOR COMPRESSOR ELECTRICAL ENGINE Fig. 1. Diagram of the experimental system Table 2. Specifications of components AAC system Components Specifications Compressor Swash-plate type 138 cc Number of Cylinders: 5 Condenser (580 x 350 x 20) mm3 Evaporator (220 x 260 x 60) mm3 Expansion valve TXV (Internally equalized with bulb 5.27 kW) Internal heat exchanger (192 x 73 x 63) mm3 Type: Brazed plate Number of plates:24 Heat transfer area: 0.6 m2 In order to provide the desired air flow temperatures in the evaporator and condenser ducts, electric heaters were used. Additionally, to provide the required air stream at the evaporator and condenser ducts, axial fans were installed. 2.1 Test Conditions The AAC system was first charged with R134a. The charge amount of R134a was adjusted based on the best COP value of the system. To obtain mass equivalence between the refrigerants, the charge amount of the R152a and R444A was determined by considering liquid density (Table 1) of the R134a. The charge amount of the R134a, R152a and R444A were 625g, 460g and 512g, respectively. First, experiments were conducted for the system with R134a. Then the similar experiments were conducted for R152a, R444A and R444A+IHX under the same conditions. Before operating the system with different refrigerants, the system was pumped down. Then, the system was loaded with the regarded refrigerant. Finally, the steady-state conditions of the system were prepared within 10 to 15 minutes. The same amount of Polyolester (POE) oil was charged in the compressor as the lubricant during the experiments. The equipment of the AAC system was fully insulated against external influences. The measuring devices were mounted, and a data acquisition system was used for transferring the measured data. Table 3 demonstrates the properties and precisions of the measuring devices. The uncertainty values for the Q , COP, W , Ex., ■z~ & s ' '-2 < b) d) Tlme [s] Fig. 3. Time domain signal: a) clean signal with noise; b) clean signal; c) intensified signal by adoption of the classical spectral subtraction; and d) Intensified signal using wavelet packet subtraction contains three components: the harmonic components related to machines, the modulation of some frequencies, and the aliasing of harmonic components. Therefore, the simulation signal of the drilling process with the sampling frequency f = 10 kHz is as follows: x(t ) = (1 + sin 5nt ) • cos(50nt + 0.2sin 10nt ) +sin200nt + n(t ), (10) where n(t) represents random noise. The signal is processed using classical spectral subtraction and WPMSS, and the three-layer wavelet a) b) X: 100.1 Y: 0.6265 X: 25.02 Y: 0.7262 100 150 Frequency [Hz] 100 150 200 Frequency [Hz] X: 25.02 Y: 0.7094 X: 100.1 Y: 0.6172 X: 235 Y: 0.1743 0 50 100 150 200 250 300 c) Frequency [Hz] Fig. 4. Spectra of the simulation signal: a) spectrum of the intensified signal by adoption of the classical spectral subtraction; b) spectrum of the intensified signal using wavelet de-noising; and c) spectrum of the intensified signal using WPMSS 0 50 200 250 300 0 50 250 300 c) 0 Feature Enhancement Method for Drilling Vibration Signals by Using Wavelet Packet Multi-band Spectral Subtraction 223 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 219-229 packet (db5) decomposition is used for frequency division processing of the simulation signal [17] to [20]. It is determined from Fig. 3 that the WPMSS in this study exceeds the classical spectral subtraction for the signals on the time domain waveform. Fig. 4 shows the spectrum of the simulation signal, intensified signal by adopting the classical spectral subtraction, wavelet de-noising and WPMSS, respectively. Fig. 4a indicates that under the influence of strong noise, in addition to the inherent frequency feature of the simulation signal 25 Hz and 100 Hz, there are also many harmonic components with large amplitudes in the signal spectrum, such as 235 Hz after classical spectral subtraction, which significantly affects the signal feature extraction. From Fig. 4b, it can be concluded that there are still some small false feature components in the simulated signal spectrum that is processed by the wavelet method. It is indicated from Fig. 4c that the WPMSS can extract the relevant features of signal components in Eq. (10) effectively, and it also has a good effect on reducing noise. In general, the SNR can be adopted to differentiate the effect exerted by signal analysis and processing [21]. The SNR is acquired as follows: SNR (yt ) = 10 x lg E s, Eh (dB), (11) where s(i) is the clean signal and y(i) is the processed signal. The SNR refers to the essential ratio of clean signal and noise, and the larger the SNR, the closer the processed signal to the clean signal. Table 1 indicates that the SNR of the simulation signal in Eq. (10) was processed with classical spectral subtraction, wavelet de-noising and WPMSS. In the simulation, it is evident that the trend of the SNR is exhibited in Table 1 by changing the clean and noise signals. This indicates that the WPMSS exceeds the classical spectral subtraction algorithm and wavelet de-noising in this study for the wide frequency band. Accordingly, the noise interference is not only eliminated, but the influence on the clean signal distortion is reduced. Table 1. SNR of the simulation signal Method Spectrum subtraction Wavelet de-noising WPMSS SNR 0.88 10.04 11.06 3 APPLICATION STUDY FOR FEATURE INTENSIFICATION WITH WPMSS IN DRILLING PROCESS 3.1 Drilling Experiments In this study, the drilling process vibration signal is selected as the research object to ascertain the effect exerted by WPMSS to intensify the features of drilling process signals. The drilling experiment is exhibited in Fig. 5. Fig. 5. Hole-drilling experiment A Kistler8793A three-way acceleration vibration sensor and a Kistler8152B acoustic emission sensor are adopted as sensors in the experiment. The sensors are all installed on the spindle of the machine tool, and the Z-axis vibration signal is selected as the object. To minimize the influence of environmental noise on the signal, the drilling process is dry cutting, and the signal acquisition equipment used is the National Instruments Corporation's NI PXle-1082 chassis. In the course of drilling, the change in contact position between the drill edge and the workpiece evidently changes the high-frequency elastic stress wave signal because of the lattice dislocation slip inside the cutting deformation area material, especially when the drill edge is cut in and out of the workpiece. Therefore, the acoustic emission (AE) signal in the experiment is adopted primarily to comparatively ascertain the monitoring signal intensification effect by adopting WPMSS. The drilling processing equipment and process parameters are listed in Table 2. First, we analyse the frequency spectrum of the drilling idling signal and the vibration signal. The spectral analysis (see Fig. 6) shows that the noise (drilling idling signal) and the noisy (drilling signal) spectrums are different in amplitude while they are occupying the same spectral band. It can be concluded that the conditions defined in the 4th paragraph of the 224 Zhou, Y. - Li, Y. - Liu, H. Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 219-229 introduction are satisfied so spectral subtraction can be applied to the drilling vibration signal. A similar analysis appears with the other researchers; they used a de-noising method by spectral subtraction to the detection of defects in ball bearings [22]. Table 2. Parameters used for experimental trials Description Values (model) Machine tool VMC-C30 Material Ti6Al4V Hole diameter 8 mm Spindle speed 600 rpm Tool feed rate 0.5 mm/s Sampling frequency 20 kHz 0.5 IP* ' " 0 20 40 60 Time [s] S136 kyi nPlw" 0 20 40 60 Time [s] J _ 0 H CL-0.2 H.............. . """""'-' .'"'"'''l -Ö. ¿^0.2 <| 0 20 40 60 Time [s] S132 0 20 40 60 Time [s] 0.2 t c^-0.2 0 20 40 60 < 0 20 40 60 Time [s] Time [s] 0 20 40 60 S137 0 mmmmmmimlimm 0 20 40 60 Time [s] Fig. 7. Time domain diagram of signal after WPD S§ -20 g -30 S c P- - Noise spectrum(idling) -Noisy spectrum(drilling) Table 3. Wavelet packet frequency range and corresponding frame 0 2000 4000 6000 8000 10000 Frequency [Hz] Fig. 6. Spectral analysis of the noise (idling) and noisy (drilling) signal 3.2 Comparison of Classical Spectrum Subtraction, Wavelet Packet and WPMSS with Drilling Signals In the 24 holes of the aforementioned drilling, the monitoring signal of one hole (No. 3 borehole) is selected randomly as the research object. It is indicated from Fig.7 that the eight different frequency band signals are attained by adopting a three-layer wavelet packet (db5) divider processing to divide the machine idle signal and the vibration signal. In addition, the appropriate frame length is determined by comparing the spectrum estimation of high frequency and low-frequency signals at different frame lengths through simulation and experiment, and the selection is based on ensuring the short-term stationarity of drilling signals in high-frequency frames and the spectrum estimation accuracy of signals in low-frequency frames. The frequency range and frame length subdivided by the wavelet packet of the vibration signal are listed in Table 3. No. Node Frequency range [Hz] Corresponding frame 1 [3,0] 0 to 1250 210=1024 2 [3,1] 1250 to 2500 29=512 3 [3,2] 2500 to 3750 29=512 4 [3,3] 3750 to 5000 28=256 5 [3,4] 5000 to 6250 28=256 6 [3,5] 6250 to 7500 27=128 7 [3,6] 7500 to 8750 27=128 8 [3,7] 8750 to 10000 27=128 3.3. Feature Intensification of Vibration Signal in the Drilling Process The clean vibration signals of No. 3 are attained via classical spectrum subtraction, wavelet packet, and WPMSS, respectively, and the time-domain waveforms of the signal before and after the intensification are exhibited in Fig. 8. It can be seen from Fig. 8b that classical spectral subtraction has no obvious de-noising effect on the signal. This is mainly attributed to the classical spectral subtraction requiring the noise to be statistically stable throughout the speech segment; however, unlike white Gaussian noise, which has a flat spectrum, the spectrum of machine idle noise is not flat. Thus, the noise signal does not affect the speech signal uniformly over the whole signal. While Fig. 8c represents the drilling signal after processing by wavelet de-noising, the weak features can be easily filtered, such as the information of the drill edge drilled in and out of the workpiece, resulting in a loss of characteristic information. Fig. 8d shows that the interference of the noise in the intensified signal by adopting WPMSS is evidently eliminated. In particular, when the drill S130 S131 S133 0 0 S135 S134 Feature Enhancement Method for Drilling Vibration Signals by Using Wavelet Packet Multi-band Spectral Subtraction 225 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 219-229 a) b) 30 40 50 60 70 Time [s] 30 40 50 60 70 Time [s] 30 40 50 60 70 Time [s] 30 40 50 60 70 d) Time [s] Fig. 8. Time domain figure of vibration signals: a) original vibration signal of the drilling process; b) intensified signal through classical spectrum subtraction; c) intensified signal through wavelet packet; and d) intensified signal through WPMSS edge is drilled in and out of the workpiece, the abrupt feature of the vibration signal is intensified, and the instantaneous vibration texture at the time of drilling into the workpiece is also remarkably evident. b) ** 20 Fig. 9. Energy distribution of 3D time-frequency domain in light of the WPD: a) wavelet packet energy distribution of noisy signal; and b) wavelet packet energy distribution of intensified signal In addition, when the vibration signal is decomposed into different frequency bands via the wavelet packet, the frequency band is confined by the energy and frequency. In this paper, the square of the WPD coefficient refers to the energy of the signal in the time domain [23]. From Fig. 9, it is evident that the energy distribution is basically the same, indicating that the method proposed in this paper has a certain inhibitory effect on noise. 3.4 The Energy Density Ratio (EDR) at the mutation before and after the signal was processed The EDR of the signal mutation and stationary period can be adopted as an evaluation index to intensify the features of the vibration signal in the drilling process. The formula for the energy density of the vibration signal can be expressed as: 1 N e=N § (12) where N is the signal sampling length and A(i) is the amplitude of the signal. The energy density ratios of the signal mutation and stationary are equated as: 0 10 20 70 70 10 20 0 10 20 c) 10 20 226 Zhou, Y. - Li, Y. - Liu, H. Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 219-229 * = ^, E„ (13) where EZ and EK are the energy densities at mutation and in the stationary period of the drilling vibration signal, respectively. For the same sampling length, the larger the K value, the more evident the mutation in the signal. In this paper, ten of the 24 holes were randomly selected, and their energy density ratios were acquired in accordance with the drill edge instantaneously and before cutting into the workpiece. Fig. 10 presents the calculation results. Table 4. EDR of mutations and stability time Number of mutation times Original signal WPMSS 1 1.652 23.421 2 1.504 22.542 3 2.189 52.441 4 1.788 30.461 5 1.315 12.901 6 1.457 17.138 7 1.495 19.655 8 1.236 17.451 9 1.877 40.322 10 1.914 39.841 These data indicate that the energy ratio of the signal mutation and stationary period that is greater than ten times larger than those before the WPMSS algorithm was adopted to address the vibration signal. 3.5 Relationship between Vibration Signal Intensification and Acoustic Emission Signal during the Post-Drilling Process in Terms of Mapping In the course of drilling, the AE signal is usually only sensitive to the stress wave generated by the microdeformation of the material, which can effectively avoid the interference of low-frequency environmental noise. Therefore, the AE is adopted in this study to ascertain the effect exerted by WPMSS on the vibration monitoring signal of the drilling process. Fig. 10 illustrates the time-domain signal and the corresponding AE signals of No. 3 drill hole, which were not processed via WPMSS. The intensified signal and the corresponding AE signals in the time domain are exhibited in Fig. 11, the intensified signal of the workpiece that has just been drilled into and the corresponding AE signals in the time domain are exhibited in Fig. 12. < 2 0 10 20 30 40 50 60 70 Time [s] Fig. 10. Time domain of untreated drilling process signals and the AE signal 0 10 20 30 40 50 60 70 Time [s] Fig. 11. Time domain of the intensified signal drilling process and the AE signal 29 29.2 29.4 29.6 29.8 30 Time [s] Fig. 12. Time domain of the stage just drilled into intensified signal and the AE signal Fig. 12 is indicative of the processed vibration signal that is basically consistent with the abrupt 8 6 4 0 8 6 4 2 0 2 & 0 Feature Enhancement Method for Drilling Vibration Signals by Using Wavelet Packet Multi-band Spectral Subtraction 227 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 219-229 change in the acoustic emission signal intensity. The results are indicative of the vibration signal of a drilling process that is perceived as 'additive noise' of the monitoring signal, and WPMSS can be adopted to eliminate the impact exerted by noise, including the machine operating signal. 4 DISCUSSION AND CONCLUSIONS (1) The theoretical calculation and simulation analysis indicate that the 'additive noise' and the monitoring signal with a broad frequency band were split into multiple sub-bands with the same frequency range by adopting WPD, and thereupon the spectral subtraction and signal reconstruction eliminated the impact exerted by the 'additive noise' on the clean signal. Accordingly, the features of the monitoring signal can be effectively intensified. (2) In line with the position change of the drill and workpiece, the monitoring signal was split into a machine idle signal and a drilling process signal. The machine idling signal acted as the 'additive noise'. On this basis, the 'spectral subtraction' was adopted to eliminate the impact exerted by environmental noise, comprising the machine idling signals on the monitoring signal of the drilling process, and the impact exerted by environmental noise, including machine idling signals, on the monitoring signal of the drilling process. (3) The experimental data indicate that the signal features of the drill as it drilled in and out of the workpiece were intensified, and the correspondence between the drilling process and the vibration signal was illuminated by adopting WPMSS, which offers more accurate identification and judgment basis for drilling processing monitoring. (4) Compared with the original spectral subtraction and wavelet transform, WPMSS can remove the signal noise more effectively, and simultaneously retain the weak features of some drilling processes. 5 ACKNOWLEDGEMENTS This work was sponsored by the National Natural Science Foundation of China (No. 51775468 and No. 51375419) and the Natural Science Founda-tion of Hunan Province (No. 2016JJ2134). 6 REFERENCES [1] Zhou, Y. (2011). New approach to Inspection of batch drilling quality based on the vibration signal mapping model. Journal of Mechanical Engineering, vol. 47, no. 7, p. 171, DOI:10.3901/JME.2011.07.171. (in Chinese) [2] Teti, R., Jemielniak, K., O'Donnell, G., Dornfeld, D. (2010). Advanced monitoring of machining operations. CIRP Annals, vol. 59, no. 2, p. 717-739, DOI:10.1016/j.cirp.2010.05.010. [3] Patra, K., Jha, A.K., Szalay, T., Ranjan, J., Monostori, L. (2017). Artificial neural network based tool condition monitoring in micro mechanical peck drilling using thrust force signals. Precision Engineering, vol. 48, p. 279-291, DOI:10.1016/j. precisioneng.2016.12.011. [4] Rimpault, X., Chatelain, J.-F., Klemberg-Sapieha, J.E., Balazinski, M. (2016). Tool wear and surface quality assessment of CFRP trimming using fractal analyses of the cutting force signals. CIRP Journal of Manufacturing Science Technology, vol. 16, p. 72-80, DOI:10.1016/j. cirpj.2016.06.003. [5] Kim, D.W., Lee, Y.S., Park, M.P., Chu, N.C. (2009). Tool life improvement by peck drilling and thrust force monitoring during deep-micro-hole drilling of steel. International Journal of Machine Tools Manufacture, vol. 49, no. 3-4, p. 246-255, DOI:10.1016/j.ijmachtools.2008.11.005. [6] Ferreiro, S., Sierra, B., Irigoien, I., Gorritxategi, E. (2011). Data mining for quality control: Burr detection in the drilling process. Computers & Industrial Engineering, vol. 60, no. 4, p. 801-810, DOI:10.1016/j.cie.2011.01.018. [7] Ferreiro, S., Sierra, B., Gorritxategi, E., Irigoien, I. (2010). Bayesian network for quality control in the drilling process. IFAC Proceedings Volumes, vol. 43, no. 4, p. 264-269, DOI:10.3182/20100701-2-PT-4011.00046. [8] Peña, B., Aramendi, G., Rivero, A., López de Lacalle, L.N. (2005). Monitoring of drilling for burr detection using spindle torque. International Journal of Machine Tools Manufacture, vol. 45, no. 14, p. 1614-1621, DOI:10.1016/j. ijmachtools.2005.02.006. [9] Ramirez, C., Poulachon, G., Rossi, F., M'Saoubi, R. (2014). Tool wear monitoring and hole surface quality during CFRP drilling. Procedia CIRP, vol. 13, p. 163-168, DOI:10.1016/j. procir.2014.04.028. [10] Xiao, W., Zi, Y., Chen, B., Li, B., He, Z. (2014). A novel approach to machining condition monitoring of deep hole boring. International Journal of Machine Tools Manufacture, vol. 77, no. 2, p. 27-33, DOI:10.1016/j.ijmachtools.2013.10.009. [11] Lee, S.H., Lee, D. (2008). In-process monitoring of drilling burr formation using acoustic emission and a wavelet-based artificial neural network. International Jour-nal of Production Research, vol. 46, no. 17, p. 4871-4888, DOI:10.1080/00207540601152040. [12] D'Addona, D.M., Conte, S., Lopes, W.N., de Aguiar, P.R., Bianchi, E.C., Teti, R. (2018). Tool condition monitoring of single-point dressing operation by digital signal processing of AE and AI. Procedia CIRP, vol. 67, p. 307-312, DOI:10.1016/j. procir.2017.12.218. [13] Hussein, R., Shaban, K.B., El-Hag, A.H. (2015). Acoustic partial discharge signal de-noising using power spectral subtraction. 228 Zhou, Y. - Li, Y. - Liu, H. Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 219-229 IEEE Conference on Electrical Insulation and Dielectric Phenomena, p. 330-333, DOI:10.1109/CEIDP.2015.7352003. [14] Kamath, S., Loizou, P. (2002). A multi-band spectral subtraction method for enhancing speech corrupted by colored noise. IEEE International Conference on Acoustics, Speech and Signal Processing, vol. 4, p. 4164, DOI:10.1109/ ICASSP.2002.5745591. [15] Paliwal, K., Wójcicki, K., Schwerin, B. (2010). Single-channel speech enhancement using spectral subtraction in the short-time modulation domain. Speech Communication, vol. 52, no. 5, p. 450-475, DOI:10.1016/j.specom.2010.02.004. [16] Deák, K., Mankovits, T., Kocsis, I. (2016). Optimal wavelet selection for the size estimation of manufacturing defects of tapered roller bearings with vibration measurement using Shannon Entropy Criteria. Strojniški vestnik - Journal of Mechanical Engineering, vol. 63, no. 1, p. 3-14, DOI:10.5545/ sv-jme.2016.3989. [17] Plaza, E.G., López, P.J.N. (2018). Analysis of cutting force signals by wavelet packet transform for surface roughness monitoring in CNC turning. Mechanical Systems and Signal Processing, vol. 98, p. 634-651, DOI:10.1016/j. ymssp.2017.05.006. [18] Smith, C., Akujuobi, C.M., Hamory, P., Kloesel, K. (2007). An approach to vibr-ation analysis using wavelets in an application of aircraft health monitoring. Mechanical Systems and Signal Processing, vol. 21, no. 3, p. 1255-1272, DOI:10.1016/j.ymssp.2006.06.008. [19] Feng, Y., Schlindwein, F.S. (2009). Normalized wavelet packets quantifiers for condition monitoring. Mechanical Systems and Signal Processing, vol. 23, no. 3, p. 712-723, DOI:10.1016/j. ymssp.2008.07.002. [20] Plaza, E.G., López, P.J.N. (2018). Application of the wavelet packet transform to vibration signals for surface roughness monitoring in CNC turning operations. Mechanical Systems and Signal Processing, vol. 98, p. 902-919, DOI:10.1016/j. ymssp.2017.05.028. [21] Pauluzzi, D.R., Beaulieu, N.C. (2002). A comparison of SNR estimation techniques for the AWGN channel. IEEE Transactions on Communications, vol. 48, no. 10, p. 16811691, DOI:10.1109/26.871393. [22] Dron, J.P., Bolaers, F., Rasolofondraibe, L. (2004). Improvement of the sensitivity of the scalar indicators (crest factor, kurtosis) using a de-noising method by spectral subtraction: Application to the detection of defects in ball bearings. Journal of Sound and Vibration, vol. 270, no. 1-2, p. 61-73, DOI:10.1016/S0022-460X (03)00483-8. [23] Han, L., Li, C.W., Guo, S.L., Su, X.W. (2015). Feature extraction method of bearing AE signal based on improved FAST-ICA and wavelet packet energy. Mechanical Systems and Signal Processing, vol. 62-63, p. 91-99, DOI:10.1016/j. ymssp.2015.03.009. Feature Enhancement Method for Drilling Vibration Signals by Using Wavelet Packet Multi-band Spectral Subtraction 229 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 230-237 © 2019 Journal of Mechanical Engineering. All rights reserved. D0l:10.5545/sv-jme.2019.2998 Original Scientific Paper Received for review: 2018-09-03 Received revised form: 2019-02-15 Accepted for publication: 2019-02-21 Research on the Surface Roughness of Hardox Steel Parts Machined with an Abrasive Waterjet Alexandru Catalin Filip - Cristin Olimpiu Morariu - Laurentiu Aurel Mihail - Gheorghe Oancea* Transilvania University of Brasov, Romania This paper presents an experimental investigation on the abrasive waterjet machining (AWJM) of Hardox steels. A full factorial plan was designed and carried out to determine how the traverse speed, the material thickness, and the material type influence the surface roughness. Two materials were machined during the experiments: Hardox 450 and Hardox 500. The experimental data were analysed using statistical methods, and a mathematical model was obtained. Additional experiments were made to validate the model. The results proved that the analysis is accurate and the mathematical model will be a useful tool in industrial environments for process planning when abrasive waterjet machining is used for the considered material. Keywords: abrasive waterjet machining, Hardox steels, surface roughness, regression model Highlights • Machining of Hardox steels is rather difficult using conventional cutting methods. • This research proved that machining Hardox steels by abrasive waterjet cutting can be an effective method for manufacturing parts made from such materials. • This paper presents the results of experimental research for revealing the correct values of the working parameters when using abrasive waterjet cutting for Hardox steels. • The research used a full factorial plan of experiments followed by a statistical analysis and a proposal of a mathematical model of the roughness variation on three parameters of the process. 0 INTRODUCTION In recent decades, one of the manufacturing processes with a significant growth rate is abrasive waterjet machining (AWJM). This was determined by some practical advantages of the method, which makes it a proven flexible manufacturing process [1] to [3]: capability of cutting almost any type of material, no thermal influence, small cutting forces, easy and fast setup, and environmentally friendly. The process is based (Fig. 1) on using highly pressurized water (1) guided through a small calibrated orifice, named "nozzle" (2), and mixed with fine abrasive particles (5) in a mixing chamber (3). The abrasive jet is then sent through a focusing tube (4), towards the surface of the material (6), which is fixed on slats (7) [1] and [2]. The high speed of the abrasive jet creates enough force to cut through the material of the part. The dynamics of the process are highly complex. Research studies have been conducted to analyse the influence of process parameters on the surface quality and the precision of parts machined by AWJM. Usually, the researchers consider the main parameters of influence [4] to [6] to be the following: the traverse speed (TS), the water pressure (WP), the abrasive flow rate (AFR) and the abrasive type, the standoff distance (STD), the material thickness (MT) and the material mechanical properties (MP). Fig. 1. Principle of abrasive waterjet cutting The TS and the WP are the most researched parameters. It is well demonstrated that for any type of part material - stainless steel [4], [6] and [7], brass [4], aluminium alloys [4] and [8] or titanium alloys [5] - the decrease of the TS improves the surface quality, but increases the manufacturing time and, consequently, the cost. However, for each type of part material, the mathematical dependence between these 230 *Corr. Author's Address: Transilvania University of Brasov, Department of Manufacturing Engineering, Bdul Eroilor nr.29, Brasov, Romania, gh.oancea@unitbv.ro Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 230-237 parameters is different and cannot be included in a single model. The surface quality of a part machined by AWJM is also one of the output parameters that was analysed in scientific works [9] to [13]. It is defined by two major characteristics: the roughness and the geometry. The effects of some process parameters (MT, TS, AFR) on the surface roughness [9] and [10] of aluminium plates were investigated. However, there was only one type of material used, and a mathematical model was not assessed. Other studies have used Taguchi [11] and the design of experiments (DOE) techniques [12] to outcome mathematical models of the selected roughness values depending on input parameters, such as TR, WP, AFR, and STD. However, both research studies stated that the models are valid only for the materials used, aluminium alloy or, respectively, specific alloyed steel. An interesting approach was made using an online vibration monitoring method [13]. The research stated that the increase of TS increases the roughness values and the vibrations signal. The mathematical models proposed can be valuable for predicting the selected roughness values according to the vibration signal collected online. However, like other studies, the models are valid only for the material used: a stainless steel. The surface geometry has also been explored [14] to [18]. Due to the dynamics of the abrasive jet of water, the cut kerf results tapered, with unparalleled walls (Fig. 2), having a so-called V-shaped taper. For the most typical shape of the kerf (Fig. 2), the dimension on the face where the cutting jet enters the material (EN_D) is higher than the one obtained at the surface where the cutting jet exits the material (EX_D). This phenomenon was mathematically modelled by both scientific research [14] and [15] and equipment manufacturers who attempted to compensate it by using a so-called "dynamic jet" [1] and [3]. However, the phenomenon remains when the equipment has a usual cutting head without dynamic tilting of the waterjet. More in-depth research was conducted to analyse the surface topology alongside the material thickness. One study [16] proposes a whole quantitative parameter that models the dependence between the input parameters (TS, MT, AFR) and the output ones (the surface roughness). The conclusions mention that further research is needed to generalize the model achieved. The topology of the cut zones was investigated in research [17] and [18], for AWJM of stainless steel parts. New methods for estimating the surface quality were proposed based on experimental analysis of the surface profile. The state of the focusing tube has been also studied [19]. The most critical parameter of influence was determined to be the surface roughness of the tube. Its increase amplifies the cavitating phenomenon, which has a negative influence on AWJM efficiency. The research previously mentioned, made in the field of AWJM, stated that the validity of the results is limited mainly to the type of the tested material and the particular conditions occurring during those tests. Consequently, if a new material is supposed to be machined by AWJ, the selection of proper technological setup to obtain a specified certain surface quality is based only on similar materials previously evaluated. If the results are not satisfactory, the trial-and-error method is recommended. This paper presents a scientific approach of the above issue when using AWJM for the material family of Hardox steels [20] and [21]. EN D ill H Thickness EX_D Taper Kerf Fig. 2. Taper kerf occurrence Manufacturing of parts made of Hardox steel alloys is usually challenging to be machined using traditional cutting processes because of their chemical composition. Such materials have high abrasion resistance and an elevated hardness. Abrasive waterjet cutting of Hardox steels has been investigated. The influence of the primary process parameters (TS, AFR, STD) was evaluated [22] to [24]. The research has demonstrated that the traverse speed is the most significant process parameter that influences the surface quality of the part. At the same time, like for other metal alloys, the best water pressure values have to be the highest possible. However, each paper considered only one value of the material thickness or one type of material; therefore, these parameters of influence were constant during the experiments. Research on the Surface Roughness of Hardox Steel Parts Machined with an Abrasive Waterjet 231 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 230-237 A comparison between plasma, laser cutting, and AWJM [25] has indicated that AWJM gives the most positive results in terms of surface quality. Based on the literature review and some of their studies [26], the authors considered the need for more comprehensive research on using AWJM for the Hardox steels family. The research was initiated within a project required by an industrial company which manufactures parts made of Hardox steels by AWJM. Their need was to determine the proper setup to obtain the specified quality parameters. A mathematical model has been developed, analysed, and validated. This model can be further used by the machine operator to select the proper value of the input parameters, so the output surface quality will result as required in practice. 1 EXPERIMENTAL SETUP AND METHOD The equipment used in the experimental activities was a Maxiem waterjet cutting machine, having a 20 HP pump and a maximum water pressure of 350 MPa. The dynamics of the machine does not allow the possibility of automatic tilting, to compensate for the V-shaped taper effect. The abrasive feeding system is a regular one, with a constant feed rate, which was carefully measured before the experiments. The main goal of the experimental research was to establish the laws of influence of the process parameters on the surface quality of parts from Hardox steels, machined with an abrasive waterjet. Three input variables were selected: A. type of material; B. traverse speed; C. part thickness. The selection of these three input variables was based on the following considerations: • water pressure directly influences the process productivity and costs; almost all studied research agrees that it is better to use the maximum values given by the equipment; • equipment used for AWJM has a constant abrasive dosing system; • previous research of the authors [26] proved that the standoff distance has a lower influence on the process, and the equipment used for experiments does not have an automatic standoff setup; • previous research of the authors proved that the main parameter is the traverse speed, as also stated by all experiments done by others. According to the above statements, the values of the constant parameters were adopted, as seen in Table 1. Table 1. Values of the constant parameters used during the tests Name of parameter Units Value Water pressure MPa 350 Abrasive flow rate g/min 350 Abrasive type mesh 80 Water nozzle diameter mm 0.279 Focusing tube diameter mm 0.838 Two types of Hardox steel were tested; they have the mechanical properties presented in Table 2 and the chemical composition in Table 3. Table 2. Mechanical properties of the materials [21] Material grade Hardox 450 Hardox 500 Hardness HBW Brinell 425 to 475 470 to 530 Yield strength [N/mm2] 1200 1300 Tensile strength [N/mm2] 1400 1550 Table 3. Chemical composition of the materials [21] Material grade Hardox 450 Hardox 500 C [%] 0.23 0.30 Si [%] 0.70 0.70 Mn [%] 1.60 1.60 P [%] 0.025 0.020 Cr [%] 1.20 1.50 The usual recommendation in such type of experiments is to plan a full factorial design [27]. Two levels were adopted for each of the three factors. The standard orthogonal matrix L8 (23) is shown in Table 4, and the values of the levels can be seen in Table 5. Those values were decided based on previous research in the field [26] and by practical reasons of the industrial company. Table 4. Matrix of the experiment, L8 (23) Run order Factor A B C 2 3 4 5 6 7 8 232 Filip, A.C. - Morariu. C.O. - Mihail, L.A. - Oancea, G. Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 230-237 Table 5. Levels of the variable parameters Code Process parameter Units Level 1 Level 2 A Material type n.a. Hardox 450 Hardox 500 B Traverse speed mm/min 30 70 C Material thickness mm 12 20 Each experiment was repeated three times for improving the statistical confidence of the data. After the trials, each part was analysed in terms of surface quality. The linear cut surface of some parts is presented in Fig. 3, showing the measuring procedure. The roughness was measured with a portable Mitutoyo SJ210 Surftest tester (Fig. 4), at the middle of the linear face of each part, approximately at 4 mm from the edge, both on the upper side and the bottom one. The values were measured on the linear region of the part because in this region the traverse speed is constant with the values according to Table 5. EN - entrance EX - exit Fig. 3. View of the cut surface with the measuring zone Fig. 4. Measuring equipment for the roughness For the evaluation of the surface roughness, the Ra parameter was used because it is by far the most common in both industry and research. The Ra parameter, being the arithmetical mean deviation of the assessed profile [28] and [29] can be considered of universal use and a very reliable parameter. The measured values of the surface roughness parameter Ra, both at the entrance and the exit sides, for the three samples of each test (P1, P2 and P3) are presented in Table 6. The arithmetic mean values were also calculated in the last column of Table 3. Table 6. Values of the measured roughness Ra on the parts Run order Side P1 P2 P3 Mean 1 entrance 2.085 2.388 2.228 2.300 exit 2.192 2.145 2.402 2.246 2 entrance 2.500 2.424 2.585 2.503 exit 2.708 2.375 2.789 2.624 3 entrance 2.414 2.283 2.431 2.376 exit 2.361 2.410 2.671 2.480 4 entrance 2.624 2.953 2.538 2.705 exit 7.095 5.208 6.060 6.121 5 entrance 2.611 2.860 2.686 2.719 exit 6.981 7.282 6.588 6.950 6 entrance 2.110 2.179 2.435 2.241 exit 2.739 2.526 2.805 2.690 7 entrance 2.193 2.103 2.231 2.175 exit 2.288 2.398 2.395 2.360 8 entrance 2.471 2.838 2.278 2.529 exit 2.737 2.794 2.833 2.788 2 ANALYSIS OF EXPERIMENTAL RESULTS As stated, the main goal of this research was to establish a mathematical model for predicting the correct value of the surface roughness depending on the input parameters during the waterjet machining of Hardox steels. To achieve the goal, a full factorial experimental plan was realized. Three input variables were considered: the traverse speed, the material thickness and the material type. For each variable, two levels were considered. The output parameter was the surface roughness. It is well known that there is a difference between the values of the surface roughness on the entrance and the exit areas of the part. This difference can be significant for certain combinations of values for the traverse speed and the part thickness. It can be stated that usually the values on the exit side are greater than the ones on the entrance. In terms of quality characteristics, the critical values of the surface roughness are on the exit side - the bottom side of the part, as they are all greater than the values on the entrance side - the upper one (see Table 5). Research on the Surface Roughness of Hardox Steel Parts Machined with an Abrasive Waterjet 233 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 230-237 The statistical analysis was done by using Minitab software [30]. An essential aspect of the statistical analysis is to quantify the amount of influence of each parameter on the process output. This is usually achieved through the analysis of means (ANOM). For the roughness values, the results are presented in Fig. 5. According to the diagrams, the traverse speed and the material thickness have a significant influence on the surface roughness compared to the material type, which can be statistically neglected. Fig. 5. Analysis of means for roughness Ra It is also essential to analyse the interaction between the input parameters. These diagrams are shown in Fig. 6. Material Speed Fig. 6. Interaction diagram for roughness Ra The amount of interaction can be evaluated by the Pareto diagram of standardized effects of the input parameters (Fig. 7). According to these diagrams, the material thickness (Factor C) and the traverse speed (Factor B) are the most significant parameters of influence, in this order of importance. The material type (Factor A) has a minimal influence. The interaction Factor B/ Factor C is important, while the interactions A/B and A/C are relatively small. Fig. 7. Pareto diagram for the roughness Ra Based on the above statements, a regression that takes into consideration the interaction B/C was analysed. The global equation is: Ra = K - 0.1311-B - 0.3095 • C + 0.01149 • B • C. (1) In Eq. (1), the factor K is a constant, and its possible values are presented in Table 7. Table 7. Values of the factor K in Eq. (1) Equation type Value of K Global 5.812 Hardox 450 5.656 Hardox 500 5.948 The statistical analysis of the global model is presented in Table 8. The R-squared (R2) value of 95.41 % and the P-value of lack-of-fit of 0.167 lead to the conclusion that the model can be used to predict the values of the process parameters. Table 8. Statistical parameters of the regression model Source DF Adj SS Adj MS P-value P-value Model 3 74.021 24.6735 138.61 0.000 Linear 2 53.738 26.8691 150.94 0.000 Speed 1 26.757 26.7569 150.31 0.000 Thickness 1 26.981 26.9812 151.57 0.000 2-Way interactions 1 20.282 20.2823 113.94 0.000 Speed* thickness 1 20.282 20.2823 113.94 0.000 Error 20 3.560 0.1780 Lack-of-fit 4 1.130 0.2826 1.86 0.167 Pure error 16 2.430 0.1519 Total 23 77.581 234 Filip, A.C. - Morariu. C.O. - Mihail, L.A. - Oancea, G. Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 230-237 3 VALIDATION OF PREDICTED MODEL The validity of the mathematical model acquired in Eq. (1) was evaluated with four experiments that were carried out in the same conditions assumed during the main experimental research. The machining was done on the same equipment, and the values of the constant parameters were those mentioned in Table 1. For each experiment, two parts were manufactured. The results are presented in Table 9. The predicted values of Ra considered in Table 9 were calculated with Eq. (1). The mean values of Ra measured after the experiments were calculated using the values of the exit side of the part since the regression model was calculated assuming the same hypothesis. The difference between the predicted and the obtained values of the roughness are presented in the last column of Table 9. These values have to be within a confidence interval, calculated based on the mean square error of the statistical data (see Table 8). The usual confidence interval (CI), corresponding to 1 - a = 95 % confidence level of each estimated effect [27], is: a CI = ±zan-r = ±0.1687. \jn (2) In Eq. (2), the factor za/2 is the 2.5 % quantile of the standardized normal distribution; a is the estimated value of standard deviation and represents the number of experimental sets of values used in the research. The analysis of data presented in Table 9 shows that all values are within the confidence interval given in Eq. (2). This fact is also proved in Fig. 8, which shows the graphical comparison between the predicted (calculated) and the measured (tested) values of the roughness Ra, for Hardox 450 steel, when the traverse speed TS = 50 mm/min. Table 9. Data of experiments made for validation Test no. 1 2 3 4 Material type 1 2 1 2 Traverse speed [mm/min] 50 50 50 50 Material thickness [mm] 20 20 12 12 Roughness Ra EX 4.391 5.281 2.837 2.529 4.996 3.792 2.234 2.248 Roughness Ra EN 2.988 2.470 2.761 2.625 2.279 2.822 2.163 2.074 Mean Ra [pm] 4.690 4.540 2.540 2.390 Predicted Ra [pm] 4.560 4.560 2.440 2.440 Variance [%] -3.0 0.4 -4.1 2.0 d 4.5 3 4 o ee 2.5 > i ^^Ra-me Ra-pre asured-1 dicted-1 / ✓ / ✓ 10 12 14 16 IS 20 Part thickness [mm] Fig. 8. Comparison between predicted and measured values of roughness Ra 4 CONCLUSIONS Experimental research has proved that special steels, like Hardox, can be machined by abrasive waterjets with real advantages on other machining processes. The main issue is to choose the proper values of the input parameters to obtain the desired ones of the output parameters. For this type of material, an orthogonal experimental plan was designed and conducted. Three input variables were selected: the traverse speed, the material thickness, and the material type. The output parameter was surface roughness. The results demonstrated that the material thickness has the most significant influence on the surface roughness. The traverse speed has the second greatest influence on the surface roughness, and the material type has minimal influence. The surface roughness increases when the material thickness increases, and when the traverse speed increases. The model of regression that predicts the influence of parameters is a nonlinear one. The global model, seen in Eq. (1), which predicts the output parameter for the two materials tested has a good confidence level, demonstrated by the experiment's validation. However, a higher level of confidence is achieved by using a specific regression for each material. It must be pointed out that this model is valid within the range of variation of the input parameters used during the experimental tests. Further research will be developed to extend the validity of the regression model for other values of the material thickness and other types of steels Research on the Surface Roughness of Hardox Steel Parts Machined with an Abrasive Waterjet 235 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 230-237 with high toughness when using the technology of abrasive waterjet machining. In the research, genetic programming algorithms will also be used [31] to optimize the AWJM parameters. 5 ACKNOWLEDGEMENTS The experiments of this research were made with the support of the bridge project PN-III-P2-2.1-BG-2016-0206, ProWaterjetTech, financed by UEFISCDI. The measurements were done on the infrastructure of the "Metrology Lab Mitutoyo" from the "Institute High Tech Products for Sustainable Development: PRO-DD", lab created with the support of Mitutoyo Romania Ltd. 6 REFERENCES [1] Olsen, J., Zeng, J. (2006). The state-of-the-art of precision abrasive waterjet cutting. Proceedings of the 8th Pacific Rim International Conference on Water Jet Technology. [2] Folkes, J. (2009). Waterjet - An innovative tool for manufacturing. Journal of Materials Processing Technology, vol. 209, no. 20, p. 6181-6189, D0l:10.1016/j. jmatprotec.2009.05.025. [3] Kong, M.C., Axinte, D.A. (2012). Capability of Advanced Abrasive Waterjet Machining and its Applications. Applied Mechanics and Materials, vol. 110-116, p. 1674-1682, D0I:10.4028/www.scientific.net/AMM.110-116.1674. [4] Akkurt, A., Kulekci, M.K., Seker, U., Ercan, F. (2004). Effect of feed rate on surface roughness in abrasive waterjet cutting applications. Journal of Materials Processing Technology, vol. 147, no. 3, p. 389-396, D0I:10.1016/j. jmatprotec.2004.01.013. [5] Hascalik, A., Caydas, U., Gurun, H. (2007). Effect of traverse speed on abrasive waterjet machining of Ti-6Al-4V alloy. Materials and Design, vol. 28, no. 6, p. 1953-1957, D0I:10.1016/j.matdes.2006.04.020. [6] Valiček, J., Hloch, S., Samardzic, I., Kuznerova, M. Zelenak, M. (2012). Influence of traverse speed on surface irregularities created by the abrasive waterjet. Metalurgija, vol. 51, no. 1, p. 43-46. [7] Loschner, P., Jarosz, K., Niestony, P. (2016). Investigation of the effect of cutting speed on surface quality in abrasive water jet cutting of 316L stainless steel. Procedia Engineering, vol. 149, p. 276-282, D0I:10.1016/j.proeng.2016.06.667. [8] Klichova, D., Klich, J., (2016). Study of the effect of material machinability on quality of surface created by abrasive water jet. Procedia Engineering, vol. 149, p. 177-182, D0I:10.1016/j. proeng.2016.06.653. [9] Begic-Hajdarevic, D., Cekic, A., Mehmedovic, M., Djelmic, A. (2015). Experimental study on surface roughness in abrasive water jet cutting. Procedia Engineering, vol. 100, p. 394-399, D0I:10.1016/j.proeng.2015.01.383. [10] Carach, J., Lehocka, D., Legutko, S., Hloch, S., Chattopadhyaya, S., Dixit, A.R. (2017). Surface roughness of graphite and aluminium alloy after hydro-abrasive machining. Advances in Manufacturing, Lecture Notes in Mechanical Engineering, p. 805-813, D0I:10.1007/978-3-319-68619-6_78. [11] Perec, A., Pude, F., Kaufeld, M., Wegener, K. (2017). Obtaining the selected surface roughness by means of mathematical model based parameter optimization in abrasive waterjet cutting. Strojniški vestnik - Journal of Mechanical Engineering, vol. 63, no. 10, p. 606-613, D0I:10.5545/sv-jme.2017.4463. [12] Ahmed, T.M., El Mesalamy, A.S., Youssef, A., El Midany, T.T. (2018). Improving surface roughness of abrasive waterjet cutting process by using statistical modelling. CIRP Journal of Manufacturing Science and Technology, vol. 22, p. 30-36, D0I:10.1016/j.cirpj.2018.03.004. [13] Hreha, P., Radvanska, A., Knapcikova, L., Krolczyk, G.M., Legutko, S., Krolczyk, J.B., Hloch, S., Monka, P. (2015). Roughness parameters calculation by means of on-line vibration monitoring emerging from AWJ interaction with material. Metrology and Measurement Systems, vol. 22, no. 2, p. 315-326, D0I:10.1515/mms-2015-0024. [14] Dumbhare, P., Dubey, S., Deshpande, Y. Andhare, A.B., Barve, P.S. (2018). Modelling and multi-objective optimization of surface roughness and kerf taper angle in abrasive water jet machining of steel. Journal of the Brazilian Society of Mechanical Sciences and Engineering, vol. 40, no. 259, D0I:10.1007/s40430-018-1186-5. [15] Wang, S., Zhang, S., Wu, Y., Yang, F. (2017). Exploring kerf cut by abrasive waterjet. The International Journal of Advanced Manufacturing Technology, vol. 93, no. 5-8, p. 2013-2020, D0I:10.1007/s00170-017-0467-y. [16] Valiček, J., Hloch, S., Kozak, D. (2009). Surface geometric parameters proposal for the advanced control of abrasive waterjet technology. The International Journal of Advanced Manufacturing Technology, vol 41, p. 323-328, D0I:10.1007/ s0170-008-1489-2. [17] Hloch, S., Valicek, J., Simkulet, V. (2009). Estimation of the smooth zone maximal depth at surfaces created by abrasive waterjet. International Journal of Surface Science and Engineering, vol. 3, no. 4, p. 347-359, D0I:10.1504/ ijsurfse.2009.027420. [18] Duspara, M., Starčevic, V., Samardžic, I. (2018). Analysis of zones created with waterjet cutting of AISI 316 L corrosion resistant steel. Tehnički vjesnik - Technical Gazette, vol. 25, no. 2, p. 616-621, D0I:10.17559/TV-20170216095042. [19] Li, D., Kang, Y., Ding, X., Wang, X., Fang, Z. (2017). Effects of nozzle inner surface roughness on the performance of self-resonating cavitating waterjets under different ambient pressures. Strojniški vestnik - Journal of Mechanical Engineering, vol. 63, no. 2, p. 92-102, D0I:10.5545/sv-jme.2016.3563. [20] Bugtacki, H., Smajdor, M. (2003). Mechanical properties of abrasion-resistant Hardox 400 steel and their welded joints. Advances in Materials Science, vol. 4, no. 2, p. 5-8. [21] SSAB AB (2018). Hardox steels, from https://www.ssab.com/ products/brands/hardox, accessed 2018-06-21. [22] Zohoor, M., Nourian ,S.H., Salehi, M. (2012). Optimum Parameters for Cutting Hard and Tough Material (Hardox Steel) by Abrasive Water Jet Cutting Process. International Journal of Advanced Design and Manufacturing Technology, vol. 5, no. 5, p. 71-76, D0I:10.1007/s00170-011-3761-0. 236 Filip, A.C. - Morariu. C.O. - Mihail, L.A. - Oancea, G. Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 230-237 [23] Servatka, M., Fabian, S. (2013). Experimental research and analysis of selected technological parameters on the roughness of steel area surface Hardox 500 with thickness 40mm cut by AWJ technology. Applied Mechanics and Materials, vol. 308, p. 13-17, DOI:10.4028/www.scientific. net/AMM.308.13. [24] Hlavac, L.M., Spadto, S., Krajcarz, D., Hlavacova, I.M. (2015). Influence traverse speed on surface quality after water-jet cutting for Hardox steel. Proceedings of the 24th International Conference on Metallurgy and Materials, p. 723-728. [25] Dahil, L., Karabulut, A., Dahil, I. (2014). Comparison of advanced cutting techniques on Hardox 500 steel material and the effect of structural properties of the material, Metalurgija, vol. 53, no. 3, p. 291-297. [26] Filip, A.C., Vasiloni, M., Mihail, L.A. (2017). Experimental research on the machinability of Hardox steel by abrasive waterjet cutting. MATEC Web of Conferences, vol. 94, art. no. 03003, DOI:10.1051/matecconf/20179403003. [27] Jiju, A. (2014). Design of Experiments for Engineers and Scientists. 2nd Edition. Elsevier, London. [28] ISO 25178-2:2012 Geometrical Product Specifications (GPS) - Surface Texture: Areal - Part 2: Terms, Definitions and Surface Texture Parameters, International Organization for Standardization, Geneva, DOI:10.3403/30154352. [29] Laperriere, L., Reinhart, G. (eds.) (2014). CIRP Encyclopedia of Production Engineering, 7th ed., Springer-Verlag, London, Heidelberg. [30] Minitab (2018). Statistical Software, from https://www. minitab.com/en-us/, accessed on 2018-10-15. [31] Kovacic, M., Brezocnik, M. (2018). Reduction of surface defects and optimization of continuous casting of 70MnVS4 steel. International Journal of Simulation Modelling, vol. 17, no. 4, p. 667-676, DOI:10.2507/IJSIMM17(4)457. Research on the Surface Roughness of Hardox Steel Parts Machined with an Abrasive Waterjet 237 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 © 2019 Journal of Mechanical Engineering. All rights reserved. D0l:10.5545/sv-jme.2018.5910 Original Scientific Paper Received for review: 2018-12-03 Received revised form: 2019-02-20 Accepted for publication: 2019-03-22 Functionality Assessment of Building a Microclimate System Utilising Solar Energy in a Cold Climate Karolis Janusevicius1* - Juozas Bielskus2 - Vytautas Martinaitis1 1 Vilnius Gediminas Technical University, Faculty of Environmental Engineering, Lithuania 2 Vilnius Gediminas Technical University, Civil Engineering Research Centre, Lithuania Bringing clean and affordable energy to the market is one of the goals of sustainable development set by the United Nations. Thermal comfort is an important aspect of efficient energy use, which plays a crucial role in ensuring health and well-being in the built environment. The majority of energy in the building sector is consumed by microclimate systems that provide thermal comfort for occupants. Design strategies, such as passive and active solar and thermal mass utilisation, reduce heat demand. When aiming to optimise thermal comfort, a reasonable combination is important in order to increase the utilisation of clean renewable energy, while conserving other resources. In this paper, a method to generate design charts is proposed to assess early design options. It aids in the selection of design parameters, based on targeted seasonal thermal comfort as a function of a complex microclimate system. In order to explore the interaction between design variables, a TRNSYS simulation model was used. An analysis of comfort conditions (based on the EN ISO 7730 method) in building spaces was performed to assess functionality. The simulation model accounted for the thermal constant in building spaces, solar utilisation and gain through glass surfaces, solar collectors and active accumulation, energy transportation, and distribution efficiency. The presented case study results showed that the lack of space heating capacity (3/4 of the calculated quantity) could be compensated for by thermal mass and a solar thermal collector without compromising thermal comfort (the percentage of people dissatisfied (PPD) was below 10 %). The highest solar fraction (36 %) was reached with the lowest fractions of space heating capacity (1/2 of the calculated quantity), due to increased demands, but this design option did not satisfy the thermal comfort conditions (PPD > 17 %). Keywords: solar thermal; early design tool; thermal comfort Highlights • In order to reduce energy consumption and environmental impacts, sustainable buildings must integrate passive and active technologies. • Integration of solar thermal utilisation to meet the energy demands of buildings is limited, due to accumulation capabilities. • The suggested method involves creating design charts to identify parameters suitable for combining passive and active technologies for solar utilisation and heat accumulation strategies to meet thermal comfort functional requirements. • The case study demonstrates how sample design charts could be generated to assess design options. • The design charts were used to identify limiting values to achieving suitable thermal comfort. 0 INTRODUCTION In September 2015, the United Nations (UN) Assembly approved the sustainable development agenda [1]. According to this agenda, member countries must take action in order to improve the sustainability of the planet. Goals up until the year 2030 were set. Multiple aspects of the 17 sustainable development goals determined by the UN were aimed at increasing the utilisation of renewable energy and improving energy efficiency [2]. It is clear that these goals could only be met by using complex measures, in terms of the most influential consumers, i.e. the built environment, due to its high environmental impact [3]. In EU households, heating and hot water alone account for 79 % (in industry 70.6 %) of total final energy consumption. Thus, it is important to focus on heating demands [4]. Approximately 84 % of heating and cooling is still generated from fossil fuels, while only 16 % is generated from renewable energy. In order to fulfil the EU's climate and energy goals, the heating and cooling sector must sharply reduce its energy consumption and decrease its use of fossil fuels [5]. Integrated design efficiency and final product quality are highly dependent on primary actions at the initial stages of any project. According to value engineering principles [6], a function could be satisfied in different ways that allow for different performances. Those aspects determine the pathway characterized by [7], which specified that alternatives should be considered at the initial stage of design. In this phase, solutions can be integrated at the lowest resource cost. In this context, the cost of considering different design alternatives becomes an important aspect of an efficient design process. When designing sustainable buildings, passive and active energy conservation solutions, in combination with renewable energy utilisation, play an important role. The primary function of a complex microclimate 238 *Corr. Author's Address: Vilnius Gediminas Technical University, Faculty of Environmental engineering, Department of Building Energetics, Sauletekio av., 11, Vilnius, Lithuania, karolis.janusevicius @vgtu.lt Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 system is to ensure comfort for the occupants. A design process sets only a portion of the parameters, but overall thermal comfort expression is typically left unaddressed. Thermal comfort is one of the aspects, along with visual, acoustic, and indoor air quality, which ensures health, productivity, and the wellbeing of people inside the building. The function of creating thermal comfort is itself composed of multiple second-order functions - ensuring sensible temperature, air velocity, and relative humidity suitable for occupant activity and clothing levels [8]. There are multiple models on how to determine thermal comfort, from the holistic person-environment systems approach to the deterministic stimulus-response models, which are mostly used in standards [9]. In cold climates, the most influential factor in terms of creating internal thermal comfort is space heating demand. Ensuring comfortable indoor environments in workspaces is one of the primary contributors to achieving better performance. Indoor environmental factors in workplaces should follow functional, as well as behavioural, requirements [10]. According to Pathak et al. [10], the functionality of the indoor environment (microclimate) can be described in terms of warmth, light, noise, space, and furniture. The most defined and parametrized among the five is thermal comfort [11] and [12]. Active and passive building microclimate subsystems influence the overall thermal comfort. One of the simplest and most effective ways of pursuing energy savings in buildings is by acting on its envelope, which represents the path for substantial heat flows [13]. A passive subsystem consists of the building envelope and the internal thermal mass. It is considered passive, because it does not consume energy to provide functionality. Active subsystems consist of heating, cooling, and ventilation systems, which consume energy to perform dedicated functions [14]. In order to most effectively utilise solar radiation, passive and active measures of heat accumulation should be used. Passive utilisation is implemented via glazing and thermal mass in the building structure. Active utilisation can be fulfilled using accumulation storage fed via thermal collectors. Other parameters, such as solar heat gain and infiltration caused by air change rate, influence the indoor environment [15]. Solar heat gains achieved through glazed building envelope elements are absorbed and accumulated in the building structure. Due to the temperature difference, heat is released from the thermal mass to the indoor space and this reduces heating energy consumption. Active microclimate systems ensure thermal comfort by supplying heat, cold air, fresh air, and light to a room. However, in this case, the focus is on space heating measures. Despite the fact that solar energy is utilised to achieve microclimate needs, total demand could not be ensured, due to an inequality between consumption and generation. Achieving equality becomes more realisable as the thermal properties of the building envelope are increased, unlocking an increase in solar utilisation and a decrease in fossil fuel use for the needs arising in a built environment. Significant daily, seasonal, and annual differences in solar radiation exist in northern countries. This variation influences the requirement for accumulation in a system that utilises solar energy for space heating and other services. It is possible to accumulate solar energy passively or actively. Passive solar utilisation is that used to heat a space without any additional energy required to run the process [16]. Solar gains collected through transparent envelope elements are accumulated in the thermal mass of a building. At night, the thermal mass at higher temperatures collects heat and reduces heating load for that period. The fact that natural thermal processes in a building may reduce the energy requirement for heating makes it possible to choose a lower design temperature for the winter [15]. Therefore, studies have already shown that the impact of different insulation-mass configurations on both heating / cooling consumption and indoor comfort varies, and often exhibits the opposite effects [17]. Physical thermal processes that take place in the building envelope may reduce the demand for heating energy [15]. Studies have shown that passive solar heating combined with energy-saving construction could reduce thermal heating requirements by up to 30 % [18]. Other studies state that the decrease varies from several percent to more than 80 % [19]. The latest research shows that by applying this combination, heating costs in the winter can be reduced by up to 58 % [20]. Passive solar heating of buildings using thermal storage in the walls has been a subject of many studies, particularly related to the Trombe wall [21] to [25]. Ordinary walls have also received research attention [13], [19] and [26] to [28]. A study [29] defined a new system combining an active-passive triple phase change material wall and a multi-surface through solar concentrators. It improved heat storage capacity of the wall interior by using the active method and could improve heat storage capacity of the phase change material (PCM) wallboards using the passive method. Other studies used PCM in cold climates to accumulate solar energy, which led to a 17.4 % annual drop in heating demand [30]. Functionality Assessment of Building a Microclimate System Utilising Solar Energy in a Cold Climate 239 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 Authors [28] analysed the simultaneous optimisation of the thickness and positions of insulation layers. The optimisation effects on the dynamic thermal characteristics of walls were also quantified. Active solar energy utilisation can be implemented by using solar collectors for water heating. Storage tank and solar collectors are the main components of any solar water heating system [31]. Heat storage is an important part of a solar collector system, which helps to manage differences in time, as well as between the use and production of heat. Study [32] analysed the contribution of solar-thermal and thermal energy storage systems to fulfilling thermal energy demand profiles. 48 case studies used 12 house types under average (smoothed) and actual (warmer) weather conditions, as well as continuous and intermittent comfort maintenance. Significant research has been conducted in regards to solar systems for space heating or hot domestic water preparation with solar collectors [30] to [33]. One study [34] presented a simplified method for optimising the key parameters of solar water heating systems based on life-cycle energy analysis, taking the phenomenon of energy mismatch into consideration. An overview of the research shows a lack of studies in the holistic assessment of systems that integrate active and passive accumulation. Specifically, those using heat storage tanks with thermal mass in building structures were lacking in the literature. In addition, active and passive solar utilisation by solar collectors and glazing were absent in the research. The goal of this article is to examine the effects of combined solutions on thermal comfort [8]. The following section presents the method that expresses the design parameters for the system, combining thermal mass and solar utilization to ensure thermal comfort. 1 METHODOLOGY In order to assess the degree of function fulfilment, the functionality (or functional efficiency) term was used. The functionality of a complex microclimate system was measured using the thermal comfort model created by Fanger [35] and adopted as standard [8]. By employing mathematical functions expressing the predicted mean vote of a statistical building occupant, the percentage of dissatisfied people (occupants) can be expressed. The Fanger thermal comfort model has a limiting value of 5 % dissatisfied people. Lower values are not possible according to the model, due to the differences in human thermal sensation. This value indicates that the function of ensuring thermal comfort for occupants is satisfied at the highest level. This paper presents an algorithm, which was used to generate design charts for the initial assessment of suitable options at the initial design phase. The algorithm expressing the workflow of design chart creation is shown in Fig. 1. Fig. 1. Workflow of design chart creation on the basis of the proposed method This tool helps designers quickly assess synergies in solar energy and heat accumulation utilisation in terms of the requirements for thermal comfort of building occupants. The execution of the provided algorithm results in design charts that can be used to predict design performance with various parameter combinations. Being able to consider a variety of design parameters and select the most appropriate one, or identify how the design should be adjusted, are the main benefits of the proposed method. 1.1 Preparation of the Design Function for Parametric Analysis In order to determine the feasibility of utilising solar energy and heat accumulation methods and to quantify the interaction with internal thermal comfort, a calculation model was developed using TRNSYS software. The model was carefully parameterized using text variables and proper output files. The model was then saved as a *.dck file. The main output parameters used to guide the design decision were 240 Janusevicius, K. - Bielskus, J. - Martinaitis, V. Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 determined on the basis of the algorithm presented in EN ISO 7730[8]. The execution of procedures was wrapped into a MATLAB function, which ran the following sequence: 1. Within MATLAB, the *.dck was copied so as not to alter the original file. 2. MATLAB searched for and replaced the values of any design variables. Simulation variables in the *.dck template file were provided as MATLAB function inputs. 3. The updated *.dck file was saved. 4. The simulation then ran in the TRNSYS environment. The updated *.dck file was executed by TRNEXE.exe using the dos(*) command in MATLAB, which was a simple command line access function. 5. After running the simulation, the results were aggregated as an output. 6. Temporary files created by the simulation were removed. This function allowed for parametric runs to be conducted to explore the field of design parameters and to create surface functions expressing the relation to thermal comfort. 1.2 Set up and Execution of Parametric Analysis The purpose of parametric analysis is to generate consistent performance data. The solution developed using this approach is of higher resolution than a single minimum or maximum point obtained using optimisation methods. The knowledge presented in a visual allows for faster understanding of the change related to a performance function in relation to input parameters. The execution of the created MATLAB function used specific inputs from the parameter grid results in the assessed thermal comfort as well as solar fraction values. The function was then executed with each parameter set in the loop, until there were no un-simulated sets left. The results of the parametric analysis were post-processed and presented as two overlaid surfaces. The MATLAB functionality was used to produce the graphs. The design charts represented the summary of simulation results. There were limited possibilities for extrapolation of design results if the explored range was limited. Due to the inconsistent behaviour of the interaction between the design parameters, it was more appropriate to expand the range of parametric analysis. Thus, extreme combinations were included in order to increase prediction accuracy. The structure of possible combinations is shown in Fig. 2. Fig. 2. Workflow of design chart creation on the basis of the proposed method The values of solar collector absorber area varied in the parametric analysis from 2.4 m2 to 7.2 m2 in increments of 1.2 m2. The thermal capacity of the building varied as a parameter. Table 1. Default values for dynamic parameters Class Simple hourly method, Cm, [J/K] Parameter value in simulation Cm, [J/K] Very light 80,000Af 1,600,000 Light 110,000Af 2,200,000 Medium 165,000Af 3,300,000 Heavy 260,000Af 5,200,000 Very heavy 370,000Af 7,400,000 The range of values used was in accordance with the international EN ISO 52016-1 [36] standard, as presented in Table 1. 1.3 Design Chart Creation In order to identify the parameters that ensured the required thermal comfort level, the results of the parametric analysis were used to create design charts. To increase the reusability of the charts, design parameters were normalised or aggregated to the derived parameters. The main comfort-influencing indicators, such as the building thermal constant, specific area of solar collectors, and heat source capacity, were selected as the main optimisation parameters. The values of these parameters were expressed using the following equations: • Building thermal constant (in accordance with the EN ISO 52016-1 standard) [36]: Cm. H ' Specific area of solar collectors: T = - (1) A ASC Af Functionality Assessment of Building a Microclimate System Utilising Solar Energy in a Cold Climate 279 (2) 241 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 The results of the parametric analysis were expressed as three-dimensional surface functions that defined the relationship between system parameters. These were used to describe system efficiency and functionality dependence on primary building and mechanical system properties. The methodology behind the chart structure is shown in Fig. 3. If these charts are viewed from the perspective of the Z-axis and altitude differences are expressed in isolines the results can be combined in monograms. These monograms can be used to predict a more efficient selection of parameters when designing a building and energy systems that utilise solar energy. Results of the parametric study were post-processed using the MATLAB software. The outcome of the study predicted the main dependencies in relation to the comfort level. Figures show the average percentage of people dissatisfied (PPD) values during the heating season. This parameter predicts the functional performance of the heating system from a comfort perspective, considering that the target average PPD value should not exceed 10 %. Fig. 3. Workflow of design chart creation on the basis of the proposed method 2 DESCRIPTION OF THE BUILDING MODEL USED FOR PARAMETRIC ANALYSIS The consistency of thermal comfort should be ensured during the entire heating season. The conditions under which the thermal comfort parameter exceeds 10 % indicate that the function of the microclimate system is not satisfied [8]. To maintain comfort for occupants, certain actions are required. The activities that are involved in creating thermal comfort (the main function) typically consume energy or other resources in order to run or create the system. The simulation model was created in the TRNSYS simulation environment and comfort level calculations were performed using a MATLAB script included in the main simulation tool. TRNSYS is a dynamic simulation tool that was developed over the course of 30 years, with flexibility in modelling systems and buildings. In TRNSYS, a model is developed within the Simulation Studio environment, which is a graphical user modelling interface. This software has a modular structure that was designed to solve complex problems in energy systems by breaking scenarios down into a series of smaller components. These components are known as "types", predefined components and algorithms that model the behaviour of common systems. The model was saved in the ASCII text format as a "deck" file (*.dck), which stores the "types" and connections between each type. The model was executed using TRNEXE, an algebraic and differential equation solver, which iteratively computes the state of the system at each time step. In this case, it was seven minutes, as it was necessary to reflect the control time constant. Because the model was stored in a text format, it could be parameterised using scripting languages, such as MATLAB. Mathematical models used to perform the simulation in TRNSYS are described in the following subsections of the paper. The macrostructure of the model can be divided into four main hydronic loops (Fig. 4): 1. Solar collector to heat storage; 2. Auxiliary heater to heat storage; 3. Heat distribution network from heat storage to heat emitting devices; 4. Heat emission system in building space. The main parameter used to measure the efficiency of the design combination was seasonal comfort, as described in the following paragraph. 242 Janusevicius, K. - Bielskus, J. - Martinaitis, V. Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 2.1 Comfort Measuring how the system functions can best be described as one's sense of comfort while being in the room. Indoor comfort indicators are defined in ISO 7730 [8]. This standard describes the predicted mean vote (PMV) and PPD indices and specifies the acceptable conditions for thermal comfort. The PMV predicts the mean value of the votes from a large group of people on the ISO thermal sensation scale (+3 = hot, 0 = neutral, -3 = cold). The PPD predicts the percentage of a large group of people likely to feel "too warm" or "too cool" [37]. PMV is a function that encompasses four environmental variables: air temperature (ta [°C]), mean radiant temperature (tmrt [°C]), relative air velocity (v [m/s]), and air humidity (i.e. vapour pressure, pa [kPa]). Activity level (i.e. metabolic rate, M [W/m2]) and clothing insulation (Icl [clo]) are also considered [38]. Thus: PMV = f (ta , tmrl, Pa , M, Id , V). (3) It is important to know the percentage of people who are dissatisfied with the environment, as it represents the most likely source of complaints. Based on experimental studies in which participants voted on thermal sensations, an empirical relationship between PMV and the predicted percentage of dissatisfied occupants (PPD) was developed as follows [38]: PPD = 100-95• exp(-0.3353PMV4 -0.219PMV2).(4) Eq. (4) indicates that even at thermal neutrality (i.e. PMV = 0), 5 % of occupants may still be dissatisfied. Comfort conditions are often tested over longer periods of time, for different types of buildings, and/ or heating, ventilation, and air conditioning (HVAC) design with the assistance of computer simulation. It is necessary to quantify the long-term comfort conditions based on an index, such that alternative designs can be compared. For these purposes, the following method is recommended in future revisions of the standards [37]. The time during which the actual PMV exceeds the comfort boundaries was considered by using a factor that is a function of the PPD. Starting with the PMV distribution on an annual basis and the relation between PMV and PPD, the following was calculated: wf = PPD„ PPD„ (5) where PPDactuai PMV is the PPD corresponding to the actual PMV and PPDpmv limit is the PPD corresponding to the PMV limit. For a characteristic period during a year, the product of the weighting factor, wf, and the time, t, were added and the result was expressed in hours: • warm period: Ewf-time, where PMV > PMVlimit, • cold period: Ewf-time, where PMV < PMVlimit. The summation of the product of "weighing factor x time" was termed "weighted time" [h]. These values were used to evaluate long term comfort conditions and acceptable weighting times of, for example, 100 h to 150 h [37]. Comfort level calculations were performed using a script run in the MATLAB environment. This process was integrated into TRNSYS parametric simulations via Type 155. Thermal environmental assessment was conducted in the B category, with PMV = ±0.5 and PPD < 10 %, where Ta = 22 ±0.5 °C. The calculation adopted the following assumptions [8]: • M = 58.15 W/m2 or 1 met • W = 0 W/m2 • Icl = 0.155 m2K/W or 1 clo • v = 0.116 m/s M, W, Icl, and v remained constant during the heating season, due to low variability and limited possibilities of predicting the variations. This was due to the uncertainty related to occupant behaviour (for M, W, I) and the complexity of air movement in the room. Other parameters, such as air temperature, radiant temperature, and relative humidity, were considered to be variables and were simulated in the model as an outcome of building and ventilation unit subsystems. The air and radiant temperatures were highly influenced by the passive and active systems. Solar radiation through glazed surfaces increased the mean radiant temperature and was accumulated in the thermal mass. Active systems influenced the air and radiant temperature through the underfloor heating and the impact on air temperature. Relative humidity varied during the heating season, due to different water vapour content in the external environment. In order to generate the variables required to assess thermal comfort, a simulation model using multiple TRNSYS components was created (see Fig. 4). The mathematical models used in these components are described in the following sections. The TESS libraries were used for model creation [39]. Functionality Assessment of Building a Microclimate System Utilising Solar Energy in a Cold Climate 243 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 2.2 Building Space Model A lumped capacitance building of Type 88 represented a single-space thermal zone with internal gains. It was selected for the study, due to the possibility of including humidity exchange as well as its computational speed, due to relatively simple heat and humidity balances. This component modelled a simple lumped capacitance single-zone structure that was subject to internal gains. It differed from other simple building models since it made no assumption about the control scheme. Furthermore, it assumed an overall U value for the entire structure. Its usefulness originates from the speed with which a building heating and/or cooling load can be added to a system simulation. While this type itself did not include solar gains, this heat balance component was added as a heat source from Type 687 [39]. This calculated the amount of solar energy and illumination transmitted through a window given only the basic information - the overall U value and solar heat gain coefficient (SHGC), as well as visual transmittance. Components of air exchange heating balance were separated into two parts: air exchange caused by infiltration and mechanical ventilation with heat recovery. Infiltration air exchange was modelled as a variable velocity when the difference in air tightness of the building at 50 Pa was equal to 0.6 h-1. This effect was included via Type 932, which used the Sherman Grimsrud Infiltration Model described in the ASHRAE Handbook of Fundamentals [40]. The model based the effects of infiltration on an effective leakage area, indoor/outdoor temperature difference, and wind speed. Type 1231 modelled low-temperature hydronic heat-distributing units, such as radiators, convectors, baseboards, and finned-tube units. These types of units supply heat through a combination of radiation and convection without fans. The mathematical model implemented in this type was based on the method presented in the ASHRAE Handbook - HVAC Systems and Equipment. [41] The combination of Type 112 and Type 667 represented fans with absolute humidity input and a sensible air-to-air heat exchanger [41]. Type 6 was an auxiliary heater that elevates the temperature of the flow stream when the temperature decreased below the set point. It operated similar to an externally controlled ON/OFF heating device. Energy is delivered to the hot storage tank via an immersed heat exchanger. The heating capacity of the auxiliary (electrical) heater varied from 15 W/m2 to 30 W/m2 and the set point temperature was set at 35 °C. imulation parameter Fig. 4. Simulation model in TRNSYS environment 244 Janusevicius, K. - Bielskus, J. - Martinaitis, V. Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 2.3 Solar Collectors Flat plate solar collectors were modelled on the basis of Type 1289, which calculated dynamic efficiency as a function of inlet temperature that could be obtained using the Hottel-Whillier equation: ij'in — Tamb ) ij'in ~ Tamb ) n =n0 - ai I ' — a2 I (6) The main properties defining solar collector specifications used for the simulation are listed below: • type, flat plate solar collectors, • optical efficiency, 80.4 %, • heat loss coefficient a1, 3.235 W/(m2 K), • heat loss coefficient a2, 0.0117 W/(m2 K2), • thermal capacity C, 5.85 kJ/(m2 K), • azimuth, 0° (south), • inclination, 65°. 2.4 Heat Storage Heat storage was modelled as a stratified storage tank with 10 nodes (Type 534). The overall tank loss coefficient was assumed to be 0.4 W/(m2K). The fluid in the storage tank interacted with the fluid in the heat exchangers (through heat transfer with the immersed heat exchangers), with the environment (through thermal losses from the top, bottom, and edges), and with the flow stream that enters and exits the storage tank. The storage tank was divided into isothermal temperature nodes. Each constant-volume node was assumed to be isothermal and to interact thermally with the nodes above and below through several mechanisms: fluid conduction between nodes and fluid movement due to destratification effects. The volume of the accumulation tank was assumed to depend on the solar collector area (ASC) and can be expressed as shown in Eq. (4). Vaccumulation 0 '07 ASC ' (7) where 0.07 is the specific accumulation tank volume per 1 m2 of the solar collector area. 3 MODEL VERIFICATION PROCEDURES For engineered systems, terminology such as "virtual prototyping" and "virtual testing" is now used in development to describe numerical simulation for the design, evaluation, and "testing" of new hardware and systems. Virtual testing is a more reasonable way to ensure proper operation and performance than building and testing a physical system. The potential legal and liability costs of failures can be staggering to a company, the environment, or the public, when building inappropriate systems without virtual testing [42]. The limitation of the proposed method is the uncertainty of predicted results. In order to control this aspect, the following validation procedures should be performed: 1. Build parametric models from verified simulation models/components. 2. Numerical errors occurring due to precision and simulation parameters should not exceed 5 % of the final result value. 3. Ensure that mass flow and energy balances are satisfied in system simulations. 4. Verify that thermal comfort parameters are in accordance with EN ISO 7730 acceptable ranges. As there was no available analytical solution or experimental data for the given system configuration, only separate component verifications could be performed. The evidence of how accurately the computational model simulates reality could thus be based on separate element verification procedures previously performed for the components used in the case studies building space [43], solar thermal collectors [44] and accumulation tank [45]. The given aspects showed that the results were sufficiently accurate for initial predictions. The next stage of the design process was to perform a more detailed simulation to measure the building performance for load sizing, energy consumption predictions, and other activities typically performed in the design process. 4 RESULTS AND DISCUSSION As an example, Fig. 4 shows the dependence of PPDavg on a specific area of solar collectors and the building thermal constant. The specific area of solar collector was the ratio between the solar collector area and the floor area. The dependence of solar fraction on the same parameters is depicted in Fig. 5. In Fig. 6, the fifth (5) surface represents the PPDavg values at full heat source capacity (30 W/m2). The fourth (4) surface represents the same result, but using 7/8 of the required heat source capacity (26.25 W/m2). The third (3) used 3/4 (22.5 W/m2), the second (2) used 5/8 (18.75 W/m2), and the first (1) used 1/2 (15 W/m2) of the required capacity. Functionality Assessment of Building a Microclimate System Utilising Solar Energy in a Cold Climate 245 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 Fig. 5. Change of PPD with varying heat source power The red plane represents the boundary where PPDavg has a constant 10 % value. Results above this boundary were considered unsatisfactory in terms of thermal comfort in the room (PPDavg > 10 %). Fig. 6. Change offsdi with varying source power In Fig. 6, the fifth surface (5) represents the solar fraction values when 1/2 of the required heat source capacity (15 W/m2) was used. The fourth surface (4) represents the same quantity at 5/8 of the required capacity (18.75 W/m2), the third (3) represents 3/4 of the required capacity (22.5 W/m2), the second (2) was 7/8 (26.25 W/m2), and the first (1) represents the full heat source capacity (30 W/m2). Fig. 7 represents the PPDavg and solar fraction when the capacity was 1/2 (15 W/m2) of the heat source required (30 W/m2) in accordance with regular sizing methods. Figs. 6 and 7 represent the same quantities as depicted in Fig. 5. The heat source capacities of 3/4 (22.5 W/m2) and 7/8 (26.25 W/m2) were used for the room space heating load. The PPDavg reached an acceptable level when the heat source capacity values were 5/8 (18.75 W/m2) and 1.0 (30 W/m2) of the required heat source capacity. These relations are shown as dashed isolines. Continuous isolines represent the solar fraction under given conditions. Fig. 7 shows that increasing the thermal constant and the specific area of solar collectors contributed to the improvement of the average PPD value for the heating season (from 14.4 % to 12.8 %). Thus, the number of unsatisfied people decreased. The PPD value varied from 12.2 % to 11.2 % when the heat source was at 5/8 of the required capacity (18.75 W/m2). The increase of specific area of solar collectors influenced the comfort level and solar fraction. Fig. 7. Change of PPD and solar fraction (fsol) when heat source power is 15 W/m2 This behaviour was quite predictable, due to the increased possibilities of solar utilisation that involve a greater solar collector area. It is important to emphasise that this combination of thermal constant and the specific area of solar collectors ensures the highest comfort level, while using the lowest amount of energy. A heat source that covered only a half of the required heating capacity (Fig. 7) could help to achieve a much higher solar fraction (0.36). 5/8 of the required capacity achieved a solar fraction equal to 0.32, but had a greater comfort level in comparison to 1/2 of the required capacity, even though the average heating season temperature was above 20 °C. In this case, the duration of discomfort (PPD > 10 %) lasted for 17 % of the heating season and, depending on design parameters, varied from 16 % to 18 %. In the second case (18.75 W/m2), the discomfort varied from 14 % to 15 % for the assessment period. Fig. 8 shows how PPDavg and fso[ depend on the same parameters when the heat source is designed to meet 3/4 of capacity. In this graph the breaking point was identified by increasing the thermal constant and the specific area of solar collectors if PPDavg was acceptable and the number of dissatisfied people dropped below 10°% (see the shaded area in Fig. 8). The highest achievable solar fraction was equal to 246 Janusevicius, K. - Bielskus, J. - Martinaitis, V. Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 0.26. Under these conditions, comfort conditions were not maintained for 14 % to 15 % of the heating season. In cases when 7/8 and full required capacity were used (30 W/m2), PPDavg reached acceptable levels under all design parameter combinations (see Fig. 9). Solar fractions of 0.26 and 0.24 were achieved, respectively. Comfort requirements were not met for 10 % to 12 % of the heating season. Fig. 8. Change of PPD and solar fraction (fso) when heat source power is 22.5 W/m2 Fig. 9. Change of PPD and solar fraction (fsol) when heat source power is 26.25 W/m2 Having expressed the PPD and solar fraction values as functions of space heating capacity and the specific area of solar collectors, a chart depicting relevant relations was produced (Fig. 10). Due to the low impact of the thermal constant, the value of PPD, restricted under given conditions (10 %), varied slightly (see the shaded area in Fig. 10). For the purposes of simplification, the dependency is shown at the thermal constant of 156 hours, which is the median value of the range examined. The isoline showing 10 % of dissatisfied people decreased in relation to the solar collector area. Thus, the lack of space heating capacity could be compensated for by an energy-efficient building that satisfies the requirements stated in the assumptions of this study. In this case, solar collectors could supply 0.14 to 0.29 of the required heat demand during the heating season. 'S 25 M Q. (O [O Q. CO -4-JL 7 ——--¡¿_ „v i O < \ \ : 0.15 035 0 2 0 25 0 3 Specific area of solar collectors H Fig. 10. PPD and solar fraction (fso) dependence on heating capacity and specific solar collector area The combination of higher thermal mass and specific area of solar collectors helped to decrease the installed heating source capacity without having a detrimental effect on thermal comfort. This study showed that the highest solar fraction was ensured when the space heating source capacity was at its lowest. Under these conditions, the thermal comfort level was not ensured (the PPD was 17 %). This combination of design parameters cannot be used, since it is not acceptable from the perspective of comfort level. In the best case presented in this paper, the comfort level was not ensured for 10 % of the heating season, while space temperature was at 22 ±2 °C. This shows that temperature was important, but it was not the only parameter since thermal discomfort can be caused by relative humidity and radiant temperature as well. From a global perspective, air velocity significantly influenced comfort level, but in this case it was assumed to be constant. A combination of these parameters led to difficulties in ensuring the comfort level throughout the heating season. Functionality Assessment of Building a Microclimate System Utilising Solar Energy in a Cold Climate 247 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 5 CONCLUSIONS In this paper, functionality was defined by the PPD, which was aggregated to reflect seasonal dissatisfaction related to thermal comfort. This approach allowed for design parameters to be determined that ensured an acceptable degree of functional satisfaction. The main conclusions of this research are as follows: 1. Separate efforts to explore passive and active heat accumulation and solar energy utilisation, as well as thermal comfort were observed in the research. The integrated approach suggested in the paper, however, has not been identified in the conducted literature analysis. 2. By using the proposed method, designers are able to create a complex solution, which takes into account multiple aspects that influence the thermal comfort level. The method establishes a connection between different design parameters (as variables) and enables the identification of synergies between different active and passive strategies in order to satisfy the function of thermal comfort. The method has the following strengths: • Generation of the design charts requires a onetime effort and the knowledge of interactions can be reused. • The method replaces intuitive rule-of-thumb methods and increases clarity and transparency in the initial design process. • The design method allows for the identification of possible trade-offs between different approaches using passive and active solar and heat accumulation utilisation techniques. • It provides a tool for sustainability-conscious designers to implement a combination of passive and active utilisation of solar energy with heat accumulation. 3. High level certainty is not required, due to the tolerated accuracy margin at the early stages of the design process. The simplifications and assumptions are considered to a certain extent to reflect the trends of systems and coupled interactions. The main sources of uncertainty in the outcome of the design method are the assumptions used to create the simulation model, as well as its validity in comparison to reality. 4. The results of the created case study showed that the estimated optimal heat source capacity was 22.5 W/m2. In relation to the specific area of solar collectors and the available building thermal constant PPDavg values, this capacity was lower than 10 %. The highest solar fraction (0.36) was achieved when the heating source capacity was at its lowest (15 W/m2), but thermal comfort requirements were not met (PPDavg is higher than 10 %), thus, this design option should be avoided. 6 FUTURE WORK The strength of the proposed method lies in the fact that it enables quick selection of building and microclimate system characteristics in the early design stages. The limitation of this method is the computational cost required to generate the design charts. The costs may increase if additional variables are included. The normalisation of design variables and the possibility of reusing the design charts creates value for the design process. Further improvements to the design chart algorithm might be achieved by including additional design variables. In order to decrease the computational cost of generating design charts, the parametric analysis could be upgraded to a combination search based on multi-objective optimisation. As the thermal comfort models evolve to include methods based on the second law of thermodynamics [46] and [47], further work may lead to design charts using exergy-based definitions [48]. 7 ACKNOWLEDGEMENTS The authors would like to express their appreciation to the Laboratory of Building Energy and Microclimate Systems (Vilnius Gediminas Technical University) for providing the TRNSYS and MATLAB software tools. 8 NOMENCLATURE Cm building thermal capacity, [J/K] H building specific heat loss, [W/K] t thermal constant of the building, [h] Af heated floor area, [m2] ASC area of solar collectors, [m2] n thermal efficiency of solar thermal collectors, [-] no optical efficiency of solar thermal collectors, [-] ax loss coefficient of solar thermal collectors, [W/(m2 K)] a2 loss coefficient of solar thermal collectors, [W/(m2K2)] Tin inlet temperature of solar collectors, [°C] Tamb ambient temperature, [°C] IT incident solar radiation, [W/m2] AT temperature, [K] ta temperature of the thermal zone, [°C] 248 Janusevicius, K. - Bielskus, J. - Martinaitis, V. Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 tmrt mean radiant temperature of thermal zone, [°C] pa vapour pressure, [Pa] M metabolic rate, [W/m2] Ic¡ clothing level of occupants, [clo] Vac volume of active accumulation tank, [m3] PPD percentage of people dissatisfied at specific time step, [%] PPDavg average seasonal percentage of people dissatisfied, [%] PMV predicted mean vote, [-] Wf weight factor for seasonal thermal comfort calculation, [-] 9 REFERENCES [1] United Nations (2015). Transforming our world: the 2030 Agenda for Sustainable Development, from https://sustainabledevelopment. un.org/post2015/ transformingourworld, accessed on 2018-10-08. [2] United Nations (2015). About the Sustainable Development Goals, from https://www.un.org/sustainabledevelopment/ sustainable-development-goals/, accessed on 2018-10-08. [3] European Commission (2018). The energy performance of buildings directive, from https://ec.europa.eu/info/sites/ info/files/epbd_factsheet_20180503_dc_v03e_final.pdf, accessed on 2018-10-08. [4] European Commission (2018). Heating and cooling, from https://ec.europa.eu/energy/en/topics/energy-efficiency/ heating-and-cooling, accessed on 2018-10-08. [5] 2010/31/EU (2016). Proposal for a Directive of the European Parliament and of the Council amending Directive 2010/31/ EU on the energy performance of buildings. European commission, Brussels, DOI:10.1007/978-1-137-54482-7_33. [6] EN 12973:2000 Value management. International Organization for standardization. Brussels. [7] Cook, R. (2007). Integrated Project Delivery a Guide. The American Institute of Architects, Chicago. [8] EN ISO 7730:2005. Ergonomics of the thermal environment - Analytical determination and interpretation of thermal comfort using calculation of the PVM and PPD indices and local thermal comfort criteria. International Organization for standardization. Geneva, D0l:10.3403/30046382. [9] Megri, A.C., El Naqa, I. (2016). Prediction of the thermal comfort indices using improved support vector machine classifiers and nonlinear kernel functions. Indoor and Built Environment, vol. 25, no. 1, p. 6-16, D0I:10.1177/1420326X14539693. [10] Pathak, P.M., Dongre, A.R., Shiwalkar, J.P. (2014). Impact of spatial, thermal and lighting parameters on the efficiency and comfort of users in Indian workspaces. Journal of Sustainable Development, vol. 7, no. 4, p. 111-123, D0I:10.5539/jsd. v7n4p111. [11] Dovjak, M., Slobonik, J., Krainer, A. (2019). Deteriorated indoor environmental quality as a collateral damage of present day extensive renovations. Strojniški vestnik -Jurnal of Mechanical Engineering, vol. 65, no. 1, p. 31-40. D0I:10.5545/sv-jime2018.5384. [12] Al horr, Y., Arif, M., Katafygiotou, M., Mazroei, A., Kaushik, A., Elsarrag, E. (2016). Impact of indoor environmental quality on occupant well-being and comfort: A review of the literature. International Journal of Sustainable Built Environment, vol. 5, no. 1, p. 1-11, D0l:10.1016/j.ijsbe.2016.03.006. [13] Leccese, F., Salvadori, G., Asdrubali, F., Gori, P. (2018). Passive thermal behaviuor of buildings: Performance of external multi-layered walls and influence of internal walls. Applied Energy, vol. 225, p. 1078-1089, D0I:10.1016/j.apenergy.2018.05.090. [14] Borel, L., Favrat, D. (2010). Thermodynamics and Energy Systems Analysis: From Energy to Exergy. Taylor and Francis Group, Boca Raton. [15] Orosa, J.A., Oliveira, A.C. (2012). A field study on building inertia and its effects on indoor thermal environment. Renewable Energy, vol. 37, no. 1, p. 89-96, D0I:10.1016/j. renene.2011.06.009. [16] Passipedia (2016). The Passive House - definition, from https://passipedia.org/basics/the_passive_house_-_ definition, accessed on 2016-06-02. [17] Stazi, F., Bonfigli, C., Tomassoni, E., Di Perna, C., Munafo, P. (2015). The effect of high thermal insulation on high thermal mass: Is the dynamic behaviour of traditional envelopes in Mediterranean climates still possible? Energy and Buildings, vol. 88, p. 367-383, D0I:10.1016/j.enbuild.2014.11.056. [18] Buker, M.S., Riffat, S.B. (2015). Building integrated solar thermal collectors - A review. Renewable and Sustainable Energy Reviews, vol. 51, p. 327-346, D0I:10.1016/j. rser.2015.06.009. [19] Aste, N., Angelotti, A., Buzzetti, M. (2009). The influence of the external walls thermal inertia on the energy performance of well insulated buildings. Energy and Buildings, vol. 41, no. 11, p. 1181-1187, D0I:10.1016/j.enbuild.2009.06.005. [20] Albayyaa, H., Hagare, D., Saha, S. (2019). Energy conservation in residential buildings by incorporating passive solar and energy efficiency design strategies and higher thermal mass. Energy and Buildings, vol. 182, p. 205-213, D0I:10.1016/j. enbuild.2018.09.036. [21] Bajc, T., Todorovic, M.N., Svorcan, J. (2014). CFD analyses for passive house with Trombe wall and impact to energy demand. Energy and Buildings, vol. 98, p. 39-44, D0I:10.1016/j. enbuild.2014.11.018. [22] Dong, J., Chen, Z., Zhang, L., Cheng, Y., Sun, S., Jie, J. (2019). Experimental investigation on the heating performance of a novel designed trombe wall. Energy, vol. 168, p. 728-736, D0I:10.1016/j.energy.2018.11.125. [23] Kontoleon, K.J., Theodosiou, Th.G. Tsikaloudaki, K.G. (2013). The influence of concrete density and conductivity on walls' thermal inertia parameters under a variety of masonry and insulation placements. Applied Energy, vol. 112, p. 325-337, D0I:10.1016/j.apenergy.2013.06.029. [24] Asan, H., Sancaktar, Y.S. (1998). Effects of Wall's thermophysical properties on time lag and decrement factor. Energy and Buildings, vol. 28, no. 2, p. 159-166, D0I:10.1016/ S0378-7788(98)00007-3. [25] Al-Sanea, S.A., Zedan, M.F. (2011). Improving thermal performance of building walls by optimizing insulation layer distribution and thickness for same thermal mass. Functionality Assessment of Building a Microclimate System Utilising Solar Energy in a Cold Climate 249 Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 238-250 Applied Energy, vol. 88, no. 9, p. 3113-3124, D0l:10.1016/j. apenergy.2011.02.036. [26] Ling, H., Chen, C., Guan, Y., Wei, S., Chen, Z., Li, N. (2014). Active heat storage characteristics of active-passive triple wall with phase change material. Solar Energy, vol. 110, p. 276-285, D0I:10.1016/j.solener.2014.09.015. [27] Guarino, F., Athienitis, A., Cellura, M., Bastien, D. (2017). PCM thermal storage design in buildings: Experimental studies and applications to solaria in cold climates. Applied Energy, vol. 185, p. 95-106, D0I:10.1016/j.apenergy.2016.10.046. [28] Alghoul, M.A., Sulaiman, M.Y., Sopian, K., Azmi, B.Z. (2009). Performance of a dual-purpose solar continuous adsorption system. Renewable Energy, vol. 34, no. 3, p. 920-927, D0I:10.1016/j.renene.2008.05.037. [29] Ampatzi, E., Knight, I., Wiltshire, R. (2013). The potential contribution of solar thermal collection and storage systems to meeting the energy requirements of North European Housing. Solar Energy, vol. 91, p. 402-421, D0I:10.1016/j. solener.2012.09.008. [30] Andersen, E., Chen, Z., Fan, J., Furbo, S., Perers, B. (2014). Investigations of intelligent solar heating systems for single family house. Energy Procedia, vol. 48, p. 1-8, D0I:10.1016/j. egypro.2014.02.002. [31] Artur, C., Neves, D., Cuamba, C.B., Leao, J.A. (2018). Comparison of two dynamic approaches to modelling solar thermal systems for domestic hot water. Sustainable Energy Technologies and Assessments, vol. 30, p. 292-303, D0I:10.1016/j.seta.2018.10.012. [32] Garnier, C., Muneer, T., Currie, J. (2018). Numerical and empirical evaluation of a novel collector storage solar water heater. Renewable Energy, vol. 126, p. 281-295, D0I:10.10.16/j.renene.2018.03.041. [33] Lugaric, L., Majdandžic, L., Škrlec, D. (2010). Countrywide positioning of domestic solar water heating systems using risk analysis and geographical information system. Strojniški vestnik - Journal of Mechanical Engineering, vol. 53, no. 1, p. 3-17. [34] Yan, C., Wang, S., Ma, Z., Shi, W. (2015). A simplified method for optimal design of solar water heating systems based on life-cycle energy analysis. Renewable Energy, vol. 74, p. 271278, D0I:10.1016/j.renene.2014.08.021. [35] Fanger, P.O. (1970) Thermal Comfort. Danish Technical Press. Copenhagen. [36] ISO 52016-1:2017. Energy needs for heating cooling, internal temperatures and sensible and latent heat loads - Part 1: Calculation procedures. International Organization for Standardization. Geneva, D0I:10.3403/30303578u. [37] Olesen, B.W., Parsons, K.C. (2002). Introduction to thermal comfort standards and to the proposed new version of EN ISO 7730. Energy and Buildings, vol. 34, no. 6, p. 537-548, D0l:10.1016/S0378-7788(02)00004-X. [38] Yang, L., Yan, H., Lam, J.C. (2014). Thermal comfort and building energy consumption implications - A review. Applied Energy, vol. 115, p. 164-173, D0I:10.1016/j. apenergy.2013.10.062. [39] Thomton, J.W., Bradley, D.E., McDowell, T.P., Blair, N.J., Duffy, M.J., LaHam, N., Naik, A.V. (2012). TESSLibs 17, vol. 6, p. 1-262. [40] ASHRAE Handbook (2001). Fundamentals. American Society of Heating, Refrigerating and Air-Conditioning Engineers (ASHRAE), Atlanta. [41] ASHRAE Handbook (2016). HVAC Systems and Equipment. American Society of Heating Refrigeration and Air Conditioning Engineers (ASHRAE), Atlanta. [42] Oberkampf, W.L, Trucano, T.G. (2002). Verification and validation in computational fluid dynamics. Progress in Aerospace Sciences, vol. 38, no. 3, p. 209-272, D0I:10.1016/ S0376-0421(02)00005-2. [43] Neymark, J., Judkoff, R., Knabe, G., Le, H.T., Durig, M., Glass, A., Zweifel, G. (2002). Applying the building energy simulation test (BESTEST) diagnostic method to verification of space conditioning equipment models used in whole-building energy simulation programs. Energy and Buildings, vol. 34, no. 9, p. 917-931, D0I:10.1016/S0378-7788(02)00072-5. [44] Plantier, C., Fraisse, G, Achard, A. (2007). Development and experimental validation of a detailed flat-plate solar collector model. Researchgate. [45] Cruickshank, C.A. (2009). Evaluation of a stratified multi-tank thermal storage for solar heating applications. PhD thesis. Queen's University. Kingston. [46] Schweiker, M., Kolarik, J., Dovjak, M., Shukuya, M. (2016). Unsteady-state human-body exergy consumption rate and its relation to subjective assessment of dynamic thermal environments. Energy and Buildings, vol. 116, p. 164-180, D0I:10.1016/j.enbuild.201.01.002. [47] Prek, M., Butala, V. (2010). Principles of exergy analysis of human heat and mass exchange with the indoor environment. International Journal of Heat and Mass Transfer, vol. 53, no. 25-26, p. 5806-5814, D0I:10.1016/j.ijheatmasstrabsfer. 2010.08.003. [48] Dovjak, M., Shukuya, M., Krainer, A. (2015). Connective thinking on building envelope - Human body exergy analysis. International Journal of Heat and Mass Transfer, vol. 90, p. 1015-1025, D0I:10.1016/j.ijheatmasstransfer.2015.07.021. 250 Janusevicius, K. - Bielskus, J. - Martinaitis, V. Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 © 2019 Journal of Mechanical Engineering. All rights reserved. D0l:10.5545/sv-jme.2018.5865 Original Scientific Paper Received for review: 2018-11-17 Received revised form: 2019-02-20 Accepted for publication: 2019-03-26 Effect of Coating Thickness on a Solid-Liquid Two-Phase Flow Centrifugal Pump under Water Medium Kaikai Luo12 - Yong Wang1* - Houlin Liu1 - Matevž Dular2 - Jie Chen1 - Zilong Zhang1 1 Jiangsu University, Research Center of Fluid Machinery Engineering and Technology, China 2 University of Ljubljana, Faculty of Mechanical Engineering, Slovenia Spraying a coating on the surface of wet parts of a solid-liquid two-phase flow centrifugal pump could effectively increase its service life. To research the effect of coating thickness on the performance of the pump, a solid-liquid two-phase flow pump with a speed of ns = 81.5 was chosen, and polyurethane coatings with different thicknesses were sprayed on the surface of the blades of the impeller. The influence of coating thickness on the performance of the pump was tested under the condition of a water medium. Combined with numerical simulation, the internal flow, pressure fluctuations and radial force of the model pump were further analysed in detail. Keeping the blade outlet angle unchanged, the head and efficiency of the pump decrease with the increase of the coating thickness, and the amplitude of the decrease increases with the increase of the coating thickness. The peak value of pressure fluctuations increases with the increase of coating thickness. The pressure value on the impeller inlet increases first and then decreases with the increase of coating thickness. Radial force of impeller increases with the increase of coating thickness, and the radial force distributes in the triangle, which is the same as the number of blades. Keywords: coating thickness, solid-liquid two-phase flow pump, numerical simulation, pressure fluctuation, radial force Highlights • A significant agreement of results under experiment and numerical simulation are captured. Additionally, the relative deviation of the head and efficiency at each flow condition is within 5 %. • The outlet pressure pulsation is concentrated in the low-frequency region, and the main peak value of the pressure fluctuation is at the blade frequency. • There is a large area of low pressure at each inlet of the impeller. As the coating thickness increases, the pressure at the impeller inlet increases first and then decreases. • With the increase of the coating thickness, the streamline in the channel becomes increasingly inhomogeneous, and the streamline in the runner changes from uniform distribution to distribution near the impeller pressure face. • The radial force of the impeller increases with the increase of the coating thickness, and the radial force of the impeller is distributed symmetrically in a triangle, which is the same as the number of blades. 0 INTRODUCTION Solid-liquid two-phase flow centrifugal pumps, indispensable equipment for solid-liquid mixed medium transportation in production and daily life, play an essential role in metallurgy, water conservation, and mining. Due to the specifics of operation conditions, the solid-liquid mixture will cause wear-and-tear damage on the wet parts of a solid-liquid two-phase flow centrifugal pump, shorten the service life, and reduce operation stability. The wear of wet parts has great influence on pump performance, such as head and efficiency reduction, vibration intensification, etc. Research studies have been made on the characteristics of erosion caused by solid-liquid flow in pumps. Noon and Kim [1] conducted three-dimensional numerical analysis for the prediction of erosion on the head and efficiency losses for lime slurry flow through centrifugal pumps and found that the temperature effect is critical for erosion damage. Lai et al. [2] employed the Eulerian-Lagrangian approach and combined it with the Grant and Tabakoff particle-wall rebound model to predict particle behaviours and erosion wear in the centrifugal pump. The results show that with the increase of the concentrations of solid particles, the frequency of impingement and rebound will increase. The most serious erosion regions are around the middle of the hub and the trailing edge of blades pressure side. Sharma [3] carried out systematic studies of the erosion wear in a centrifugal slurry pump. The results show that the blade erosion rate is 10 times greater than the hub erosion rate. Compared with other parts of the impeller, the erosion rate near the leading edge is more severe. Bross and Addie [4] developed a simple model and predicted the influence of the different parameters on the wear behaviour of an impeller nose in the centrifugal slurry pump. The results show that the wear of the impeller and the liner will decrease with the increase of pressure coefficient and clearing vane height. *Corr. Author's Address: Jiangsu University, Research Center of Fluid Machinery Engineering and Technology, Zhenjiang, China, wylq@ujs.edu.cn 251 Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 In addition, some scholars researched the performance of solid-liquid flow in pumps. Wu et al. [5] researched three-dimensional unsteady flow characteristics in the centrifugal slurry pump for solid-liquid flow. The results show the fluctuation of the pressure and flow velocity, particularly near the volute tongue. The fluctuations of the head and radial forces are 8.1 % and 85.7 %, respectively. Liquids are more influenced by the clocking effect than solids are. Shi et al. [6] designed a new experimental facility without agitation for an internal solid-liquid two-phase flow test in a centrifugal pump with PIV. Ning and Wang [7] studied the solid-liquid mixed flow hydraulic characteristics and turbulent kinetic energy distribution in centrifugal pumps. Wang and Qian [8] found that the head and shaft power decrease as the silt concentration and silt particle size increase for a double-section pump with experimental methods. Mishra and Ein-Mozafari [9] analysed the performance of the Maxblend impeller in solidliquid operation using the Eulerian-Eulerian method, standard k- s turbulence model and sliding mesh technique. Cando et al. [10] adopted the filter-based RANS (k - s) and Eulerian-Lagrangian method to analyse unsteady liquid-solid flow around a step in a rectangular channel, and their results revealed that the vorticity described by the Q-criterion promotes the particle motion at the wake area. Zhao et al. [11] researched the effect of the diameter size of particles and concentration on external performance, pressure and turbulent kinetic energy distribution. Generally, coating material with better wear resistance, including epoxy resin mortar coatings [12] to [14], composite nylon coatings [15], polyurethane coatings [16] and other polymer coatings [17], sprayed on the surface of wet parts, is used to deal with the problem of the wear and tear of solid-liquid two-phase flow pumps. After spraying the coating, parameters such as the thickness and wrapping angle of pump impeller are changed, which has a significant influence on pump performance. It is necessary to focus on the effect of coating thickness on pump performance. In this work, a solid-liquid two-phase flow pump with a specific speed ns = 81.5 was chosen as the research subject. The surface of the impeller of the pump was sprayed with polyurethane coatings of different thickness to research the effect of coating thickness on the performance of the model pump. 1 RESEARCH SUBJECT, SCHEME AND TEST SYSTEM 1.1 Research Subject A two-phase flow centrifugal pump with cylindrical blades was chosen as the research subject, whose nominal flow rate Qn = 20 m3/h, head H = 22 m and efficiency n = 48 %. The main parameters of the model pump are summarized in Table 1, and its photo is shown in Fig. 1. Table 1. Main parameters of the model pump Specific speeds ns = 81.5 Nominal flow rate Qn = 20 m3/h Rotational speed n = 2900 r/min Impeller inlet diameter Di = 65 mm Impeller outlet diameter D2 = 160 mm Outlet impeller width b = 7 mm Thickness of blade S = 6 mm Number of blades z = 3 Base diameters of volute D3 = 182 mm Inlet width of volute b3 = 26 mm Diffuser outlet diameter D4 = 50 mm Fig. 1. Photo of the model pump with a cylindrical blade 1.2 Test Rig and Measuring Point Distribution The test system included the model pump, pressure fluctuation testing system, electromagnetic flowmeter and power supply control system, and other features. The pump performance was tested in the test rig schematized in Fig. 2a. The water was pumped from and returned to a 10 m3 reservoir. The shaft torque and rotational speed were monitored by a torque and speed sensor with an error rate of ±0.10 %. Static pressure values were measured at the inlet and outlet of the pump by a differential pressure transfer, and the uncertainty was within ±0.10 %. The flow rate was 252 Luo, K. - Wang, Y. - Liu, H. - Dular, M. - Chen, J. - Zhang, Z. Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 Fig. 2. Test system and measuring point distribution; a) test system of the model pump: 1. outlet valve, 2. inlet valve, 3. electromagnetic flowmeter, 4. model pump, 5. performance test system, 6. pressure fluctuation testing system, 7. power control system; b) distribution of measuring point on static pressure and pressure fluctuation measured by using a magnetic flow meter with an uncertainty of ±0.14 %. The measurement accuracy of pump efficiency was quantified as ±0.30 %. The inlet and outlet pressure measurement points were set at twice the diameter of the inlet and outlet flanges. The pressure fluctuation measurement point was set at 4 times the outlet diameter of the outlet flange. The distribution of measurement points is shown in Fig. 2b. 1.3 Scheme Non-dimensional thickness of coating thickness was dealt with to analyse the effect caused by coating thickness expediently, whose ratio was defined as K, r, Fig. 3. Coating sprayed on impeller; a) coating thickness ratio K0; b) coating thickness ratio K1, c) coating thickness ratio K2; d) coating thickness ratio K3 Effect of Coating Thickness on a Solid-Liquid Two-Phase Flow Centrifugal Pump under Water Medium 253 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 -S. ! S where Ki is the ratio of coating thickness. Si represents coating thickness, (i = 0, 1, 2, 3), and S is the thickness of the blade. The polyurethane material with better comprehensive properties was sprayed on the pressure surface of the impeller. The ratio of coating thickness is K0 = 0 (S0 = 0 mm), K1 = 1/6 (S1 = 1 mm), K2 = 1/3 (S2 = 2 mm) and K3 = 1/2 (S3 = 3 mm). The pictures of coating sprayed on impeller is shown in Fig. 3. 2 COMPUTATION MODEL AND NUMERICAL SIMULATION 2.1 Computation Model The computation model, including the impeller, volute, inlet extension and outlet extension, was built with Creo 3.0 commercial software. The purpose of adding an extension of the inlet and outlet is to reduce velocity gradient and improve simulation accuracy. The computational domain of model pump with K0 is shown in Fig. 4. 2.2 Computational Domain Discretization and Mesh Sensitivity Analysis The fluid domain of the model pump was discretized by the commercial software ANSYS ICEM. To guarantee fluid development enough, the length of the inlet extension is 5 times the inlet diameter of the centrifugal pump, and the outlet extension length is 10 times longer than the outlet diameter of the volute. Unstructured tetrahedral cells were used to discretize the computational domain. Computational model grids of the model pump are shown in Fig. 5. A detailed mesh sensitivity analysis of the computational domain was performed to eliminate the influence of the mesh factor. The same topological structure of the computational domain was employed for the research model. Seven schemes of grids were obtained at the same premise of mesh quantity. The numbers of grids for seven cases were 975,000, 1,203,000, 1,415,000, 1,694,000, 1,981,000, 2,228,000 and 2,602,800. Fig. 4. Computational domain of model pump with K0; 1 inlet extension; 2 volute 3 impeller; 4 outlet extension Fig. 6. Mesh sensitivity analysis The head value was taken as the criterion to evaluate the mesh sensitivity. As can be observed in Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 Fig. 6, the head is gradual decline with the number of girds increasing. Relative deviation of head value between case 6 and case 7 was within ±1 %, which indicates that it would be little impact on numerical results to increase the number of grids from those of case 6. Taking the computer resources in account and keeping the cost under control, the girds numbers of case 6 were chosen. 2.3 Boundary Condition and Solution Method The velocity inlet and pressure outlet boundary conditions of the model pump were loaded. The hydraulic diameter D of the inlet of the model pump is 65 mm and the value of turbulence intensity I was calculated directly [18]. Sliding mesh technology was employed to deal with the information transfer between the rotating impeller and the stationary volute. The walls of the model pump were defined as no-slip condition, and the standard wall function was used to solve the low Reynolds number flow near the wall region. Commercial software Fluent based on the finite volume method was amplied to simulate unsteady flow field of model pump. ReNormalisation Group (RNG) k-s turbulence model and the SIMPLEC algorithm were chosen to solve the turbulence function and pressure velocity coupling, respectively. In order to obtain detailed resolution of unsteady flow results of the solid-liquid two-phase flow centrifugal pump, the time step of the unsteady calculation was set as At = 5.747*10-5, which indicate that each impeller revolution will be calculated in a time sequence of 360 times steps corresponding to 1° of the impeller rotation speed. The numerical residual convergence criterion was set as 10-5 so as to ensure the result to be converged. To obtain a rapid convergence process, steady simulation was first carried out in advance, which was set as the initial condition for transient calculation. 3 RESULTS AND DISCUSSIONS 3.1 Performance Characteristics Analysis Comparison of the performance curves of four schemes are shown in Fig. 7. With the increase of flow rate, the heads value of experiment and numerical results for four schemes decrease, and the efficiency value of four schemes increases first and then decreases. The experiment results show that head and efficiency decrease with the increase of the coating thickness, but their rate of decline increases. Compared with the model pump at nominal flow rate, the head decreases by 3.74 %, 10.47 % and 24.24 % respectively when the coating thickness coefficients are Kj, K2 and K3, and the efficiency decreases by 3.57 %, 6.22 % and 14.88 % respectively. It indicates that the coating thickness has great influence on the performance of the model pump. It is mainly because the coating increases the actual thickness of the blade and reduces the flow area of the impeller passage. The performance of the pump is changed. However, the effect of the coating thickness on the performance of the model pump is smaller under low flow rate. ----_ __" VS. • s N * VV 's S.vVVv's. - >—Ov • \ N. SN.' \ • N S N^T^y^ . -Kq (Experiment) — -K2 (Experiment) ^ ---Kq (Simulation) ---K2 (Simulation) -K\ (Experiment) -Kt, (Experiment) ---K\ (Simulation) ---Kt, (Simulation) a) b) QIQn Fig. 7. Mesh sensitivity analysis; a) head curves; b) efficiency curves A significant agreement of results between experiment and simulation is captured. Under the same operating condition, the head and efficiency decrease gradually with the increase of the coating thickness. The larger the coating thickness, the greater the decline amplitude. Relative deviations of the head and efficiency at each flow condition are within 5 %. Effect of Coating Thickness on a Solid-Liquid Two-Phase Flow Centrifugal Pump under Water Medium 255 Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 3.2 Pressure Pulsation Analysis The result of the transient pressure value at the pump outlet was obtained. To make the comparion accurately, the pressure coefficient is introduced: where p is the transient static pressure of pump outlet. p is the average static pressure for a rotation period of impeller and u2 is the circumferential velocity of the impeller. The experimental results of the pressure fluctuation frequency spectrum of four schemes under different coating thickness are shown in Fig. 8, which 0 2 4 6 8 10 0 2 4 6 8 10 c) f,fn d) f,f* Fig. 8. Pressure fluctuation frequency spectrum; a) K0, b) K1, c) K2, and d) K3 were obtained by fast Fourier transform (FFT). f is the frequency obtained by FFT. f is the rotation frequency of the pump shaft. It can be seen from Fig. 8 that the trend of the pressure fluctuation of the outlet with different coating thickness is basically the same under different operating conditions. The pressure pulsation is found in the low frequency region. Due to the rotor-stator interaction between the impeller and the volute, the main peak value of the pressure fluctuation is at the blade passing frequency. When the coating thickness is the same, the pressure fluctuation coefficient under off-design condition is larger than that under the design condition. It indicates that the flow stability of the off-design condition is worse than that of design condition. The amplitude of the pressure pulsation frequency spectrum of the numerical simulation after FFT is Fig. 9. Peak value of pressure fluctuation between experiment and numerical 256 Luo, K. - Wang, Y. - Liu, H. - Dular, M. - Chen, J. - Zhang, Z. Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 extracted and compared with the experimental value as shown in Fig. 9. It can be seen from Fig. 9 that the variation trend of the numerical calculation and experimental results under different operation conditions is basically the same with the same coating thickness. The peak value of pressure fluctuation at nominal condition is smaller than that at off-design condition. The time-domains of the numerical simulation pressure fluctuation under different operating conditions in one impeller rotation cycle with diffident coating thickness are shown in Fig. 10. It can be seen from Fig. 10 that the variation trend of pressure fluctuation time-domain is identical under different operation conditions with the same coating thickness. Three peaks and three troughs could be captured in a rotation cycle, which is the same as the number of impeller blades. It indicates that rotorstator interaction is the main reason of the pressure pulsation. When the operation condition remains the same, the peak value of the pressure fluctuation of scheme K1 is the smallest. With the increase of the coating thickness, the peak value of the pressure fluctuation increases gradually. It indicates that the a) bj Anglen Cj Fig. 10. Time-domain of numerical simulation pressure fluctuation; a) 0.8 Qn, b) 1.0 Qn, and c) 1.2 Qn al) bi) clj c2j ^ c3j - c4j Fig. 11. Pressure distributions at the middle section of impeller; al) K0, 0.8 Qn; a2) K1, 0.8 Qn; a3) K2, 0.8 Qn; a4) K3, 0.8 Qn; bl)Ko, 1.0Qw b2)Ki, 1.0Qn, b3)K2, 1.0Qn; b4)K3 , 1.0Qn; cl)Ko, 1.2Qn; c2)K , 1.2Qn; c3)K2, 1.2Qn; c4)K3 , 1.2Qn Effect of Coating Thickness on a Solid-Liquid Two-Phase Flow Centrifugal Pump under Water Medium 257 Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 coating thickness of blade will increase the instability of pump outlet flow. The thickness of the blade will be increased after the coating treatment, which results in the enhanced rotor-stator interaction between the impeller and the tongue of volute. 3.3 Internal Flow Analysis Fig. 11 shows the pressure distributions at the middle section of the impeller under the conditions of 0.8 Qn, 1.0 Qn and 1.2 Qn, respectively. As can be easily seen in Fig. 11. The variation trend of the static pressure is basically the same under different operation conditions with the same coating thickness. Low pressure area at the inlet of the impeller is captured. With the increase of the flow rate, the range of the low pressure area increases gradually. With the increase of the coating thickness, the pressure at the inlet of the impeller increases first and then decreases under high flow rate and design conditions, which is related to the change of the relative flow angle of fluid. It is clarified that the coating thickness on impeller should be within a certain range, and the static pressure at the inlet of the impeller can be increased relatively within this range. Based on the above analysis, it is considered that the static pressure at the inlet of the impeller can be increased when the coating thickness coefficient is K~K2 under the condition of clear water. Fig. 12 mainly focus on the relative velocity and streamline diagram at the middle section of the impeller under the conditions of 0.8 Qn, 1.0 Qn and 1.2 Qn respectively. It can be seen from Fig. 12 that under the same operation condition, different sizes of low speed zones in the impeller runner with different coating thickness are caught, and the low speed zones with uneven distribution will lead to the instability of flow in the runner. With the increase of the coating thickness, the area of the low-speed area in the impeller runner decreases first and then increases. It is shown that proper coating thickness on blades can reduce the ] Pressure [MPa] 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 1.3 alj b1) cl) ^ v c2j c3) c4) Fig. 12. Velocity and streamline diagram at middle section of Impeller;; a1) K0, 0.8 Qn; a2) Kx, 0.8 Qn; a3) K2, 0.8 Qn; a4) K3 , 0.8 Qn b1) K0 , 1.0 Qw b2) K , 1.0 Qn b3) K2 , 1.0 Qn, b4) K3 , 1.0 Qn, Cl) K0 , 1.2 Qn c2) Kj , 1.2 Qn, c3) K2 , 1.2 Qn, c4) K3 , 1.2 Qn 258 Luo, K. - Wang, Y. - Liu, H. - Dular, M. - Chen, J. - Zhang, Z. Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 instability in the runner, which is consistent with the analysis of pressure fluctuation. The low speed zone decreases gradually with the increase of the flow rate under the same thickness. As can be shown in Fig. 12 approximate circular vortex could be found near the impeller inlet of four scheme in each operation condition, and the vortex velocity is small. It is since that the model pump impeller does not have a front cover, and a vortex is formed by the impeller rotating. Higher liquid velocity outside the vortex at the impeller inlet results in lower pressure here, which is consistent with the static pressure analysis above. Under the same operation condition, the circular vortex at the inlet decreases with the increase of the coating thickness, which leads to the increase of the inlet velocity of the impeller. The streamline distribution in the channel gradually changes from uniform to near the impeller working face, and the streamline near the low pressure region has similar vortex. When the coating thickness is the same, the range of the low pressure zone under the design condition is smaller than that under non-designed condition, and the streamline distribution is more uniform, which indicates that the coating treatment does not change the inherent property of the model pump. 3.4 Radial Force Analysis of Impeller Fig. 13 is the radial force distribution of the pump impeller for four coatingschemes. It can be seen from Fig. 13 that the radial force of the impeller varies periodically with time under various working conditions and distributes symmetrically in a triangle, which is the same as the number of blades. Under the same operation condition, the radial force increases with the increase of the coating thickness. It is because that the flow area decreases after the coating treatment, which result in the change of the flow stability in the passage. The radial force of impeller is minimum at the designed condition under the same coating thickness. Fig. 14 is the histogram of the radial force amplitude of the pump impeller with four coating thickness models. Compared with the coating thickness coefficient K0, the radial force of the impeller increases by 14.28 %, 22.06 % and 32.31 % when the coating thickness coefficients are K1, K2 and K3, respectively. -k0 -k\ -k2- -k3 120 100 80 60 150 / /......... 40 /7 20 0 180 20 \ 40 V; 60 80 210 100 120 a) 30 330 b) c) Fig. 13. Radial force of Impeller; a) 0.8 Qn; b) 1.0 Qn; c) 1.2 Qn Fig. 14. Comparison of radial force for four schemes Effect of Coating Thickness on a Solid-Liquid Two-Phase Flow Centrifugal Pump under Water Medium 259 Strajniski vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 4 CONCLUSIONS In this paper, polyurethane coating with good wear resistance was sprayed on the pressure surface of the impeller blade. The performance, pressure fluctuation, internal flow and radial force of four kinds of coating schemes under water condition were studied by experiment and numerical simulation. 1) The experimental results show that the head and efficiency decrease with the increase of the coating thickness, and the decline amplitude is proportional to the coating thickness 2) The outlet pressure pulsation is concentrated in the low frequency region, and the main peak value of pressure fluctuation is at the blade frequency. Under the design condition, the amplitude and peak value of pressure fluctuation are the smallest when the coating thickness coefficient is K1. The peak value of pressure fluctuation increases gradually with the coating thickness increasing. 3) There is a large area of low pressure at each impeller inlet. As the coating thickness increases, the pressure at the impeller inlet increases first and then decreases. 4) There are different sizes of low velocity zones in impeller flow channel under each coating thickness. With the increase of the coating thickness, the streamline in the channel becomes more and more inhomogeneous, and the streamline distribution in the impeller changes from uniform to near the impeller pressure surface. 5) The radial force of the impeller increases with the increase of coating thickness, and the radial force of impeller distributes symmetrically in triangle, which is the same as the number of blades. 5 ACKNOWLEDGEMENTS The authors would like to thank the support by the National Key Research and Development Program of China (2017YFC0804107), the National Natural Science Foundation of China (51779106, 51609164), Ministry of Education, Xihua University (szjj2016-068), Jiangsu top six talent summit project (GDZB-017), the Natural Science Foundation of Jiangsu Province (Grant No. BK20160574). 6 NOMENCLATURES Qn nominal flow rate, [m3/h] H head, [m] n rotational speed, [rpm] ns specific speeds, ns = 3.65nQ1/2/H3/4 n efficiency u2 impeller circumferential velocity, [m/s] D impeller inlet diameter, [mm] D2 impeller outlet diameter, [mm] D3 inlet diameters of volute, [mm] D4 diffuser outlet diameter, [mm] b3 inlet width of volute, [mm] z blade number of impeller p average static pressure in a rotation period of impeller, [Pa] 8 thickness of blade, [mm] ¿j coating thickness p density, [kg/m3] b outlet impeller width, [mm] f rotation frequency of pump shaft, [Hz] K ratio of coating thickness Cp time-dependent nondimensional pressure coefficient 7 REFERENCES [1] Noon, A.A., Kim, M.-H. (2016). Erosion wear on centrifugal pump casing due to slurry flow. Wear, vol. 364-365, p. 103111, D0l:10.1016/j.wear.2016.07.005. [2] Lai, F., Zhu, X., Xu, X., Li, G. (2018). Erosion wear and performance simulation of centrifugal pump for solid-liquid flow. ASME Power Conference collocated with the ASME 12th International Conference on Energy Sustainabillty and the ASME Nuclear Forum, p. V002T10A001-V002T10A001, D0I:10.1115/P0WER2018-7151. [3] Sharma, A.K., (2008). Numerical Study of Erosion Wear on Centrifugal Slurry Pump. MSc thesis, Thapar University, Patiala. [4] Bross, S., Addie, G. (2002). Prediction of impeller nose wear behaviour in centrifugal slurry pumps. Experimental Thermal and Fluid Science, vol. 26, no. 6-7, p. 841-849, D0I:10.1016/ S0894-1777(02)00174-7. [5] Wu, B., Wang, X.-L., Liu, H., Xu, H.-L. (2015). Numerical simulation and analysis of solid-liquid two-phase three-dimensional unsteady flow in centrifugal slurry pump. Journal of Central South University, vol. 22, no. 8, p. 3008-3016, D0I:10.1007/s11771-015-2837-7. [6] Shi, B.C., Wei, J., Zhang, Y. (2017). A novel experimental facility for measuring internal flow of Solid-liquid two-phase flow in a centrifugal pump by PIV. International Journal of Multiphase Flow, vol. 89, p. 266-276, D0I:10.1016/j. ijmultiphaseflow.2016. 11.002. [7] Ning, C., Wang, Y. (2016). Performance analysis on solidliquid mixed flow in a centrifugal pump. IOP Conference Series: Materials Science and Engineering, vol. 129, no. 1, p. 012062, D0I:10.1088/1757-899X/129/1/012062. [8] Wang, Z., Qian, Z. (2017). Effects of concentration and size of silt particles on the performance of a double-suction centrifugal pump. Energy, vol. 123, p. 36-46, D0I:10.1016/j. energy.2017.01.142. 260 Luo, K. - Wang, Y. - Liu, H. - Dular, M. - Chen, J. - Zhang, Z. Strojniski vestnik - Journal of Mechanical Engineering 65(2019)4, 251-261 [9] Mishra, P., Ein-Mozaffari, F. (2017). Using computational fluid dynamics to analyze the performance of the Maxblend impeller in solid-liquid mixing operations. International Journal of Multiphase Flow, vol. 91, p. 194-207, D0l:10.1016/j. ijmultiphaseflow.2017.01.009. [10] Cando, E., Yu, A., Zhu, L., Liu, J., Lu, L., Hidalgo, V., Luo, X.W. (2017). Unsteady numerical analysis of the liquid-solid two-phase flow around a step using Eulerian-Lagrangian and the filter-based RANS method. Journal of Mechanical Science and Technology, vol. 31, no. 6, p. 2781-2790, D0I:10.1007/ s12206-017-0521-6. [11] Zhao, W.G., Zhao, G.S. (2018). Numerical investigation on the transient characteristics of sediment-laden two-phase flow in a centrifugal pump. Journal of Mechanical Science and Technology, vol. 32, no. 1, p. 167-176, D0I:10.1007/s12206-017-1218-6. [12] Golru, S.S., Attar, M.M., Ramezanzadeh, B. (2014). Studying the influence of nano-Al2O3, particles on the corrosion performance and hydrolytic degradation resistance of an epoxy/polyamide coating on AA-1050. Progress in Organic Coatings, vol. 77, no. 9, p. 1391-1399, D0I:10.1016/j. porgcoat.2014.04.017. [13] Esteves, M., Ramalho, A., Ferreira, J.A.M., Nobre, J.P. (2013). Tribological and mechanical behaviour of epoxy/nanoclay composites. Tribology Letters, vol. 52, no. 1, p. 1-10, D0I:10.1007/s11249-013-0174-2. Effect of Coating Thickness on a Solid-Liquid Two-Phase Flow Centrifugal Pump under Water Medium 261 [14] Chairman, C.A., Kumaresh Babu, S.P. (2013). Mechanical and abrasive wear behavior of glass and basalt fabric-reinforced epoxy composites. Journal of Applied Polymer Science, vol. 130, no. 1, p. 120-130, D0l:10.1002/app.39154. [15] Jin, J., Rafiq, R., Gill, Y.Q., Song, M. (2013). Preparation and characterization of high performance of graphene/nylon nanocomposites. European Polymer Journal, vol. 49, no. 9, p. 2617-2626, D0I:10.1016/j.eurpolymj.2013.06.004. [16] Kim, B.K., Lee, J.C. (1996). Waterborne polyurethanes and their properties. Journal of Polymer Science Part A: Polymer Chemistry, vol. 34, no. 6, p. 1095-1104, D0I:10.1002/(SICI)1099-0518(19960430)34:6<1095::AID-P0LA19>3.0.C0;2-2. [17] Fu, J., Chen, T., Wang, M., Yang, N., Li, S., Wang, Y., Liu, X. (2013). Acid and alkaline dual stimuli-responsive mechanized hollow mesoporous silica nanoparticles as smart nanocontainers for intelligent anticorrosion coatings. ACS Nano, vol. 7, no. 12, p. 11397-11408, D0I:10.1021/ nn4053233. [18] Zhang, N., Yang, M., Gao, B., Li, Z., Ni, D. (2016). Investigation of rotor-stator interaction and flow unsteadiness in a low specific speed centrifugal pump. Strojniški vestnik -Journal of Mechanical Engineering, vol. 62, no. 1, p. 21-33, D0I:10.5545/sv-jme.2015.2859. Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4 Vsebina Vsebina Strojniški vestnik - Journal of Mechanical Engineering letnik 65, (2019), številka 4 Ljubljana, april 2019 ISSN 0039-2480 Izhaja mesečno Razširjeni povzetki (extended abstracts) Andrej Jeromen, Edvard Govekar: Nelinearni model laserskega tvorjenja kapljic iz kovinske žice na osnovi dinamičnega ravnovesja dušenega vzmetnega nihala SI 27 Mehmet Direk, Mehmet Selçuk Mert, Eren Soylu, Fikret Yuksel: Eksperimentalna raziskava avtomobilske klimatske naprave s hladilnima sredstvoma R444a in R152a kot alternativo za R134a SI 28 Youhang Zhou, Yong Li, Hanjiang Liu: Metoda za ojačenje značilk v signalu vibracij pri vrtanju z večpasovnim spektralnim odštevanjem po valčni paketni dekompoziciji SI 29 Alexandra Catalin Filip, Cristin Olimpiu Morariu, Laurentiu Aurel Mihail, Gheorghe Oancea: Raziskava površinske hrapavosti delov iz jekel Hardox, obdelanih z abrazivnim vodnim curkom SI 30 Karolis Januševičius, Juozas Bielskus, Vytautas Martinaitis: Ocena funkcionalnosti sistema za ustvarjanje mikroklime v stavbah v hladnih podnebjih s sončno energijo SI 31 Kaikai Luo, Yong Wang, Houlin Liu, Matevž Dular, Jie Chen, Zilong Zhang: Vpliv debeline prevleke na tok tekočine in trdnih delcev v centrifugalni črpalki SI 32 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, SI 27 © 2019 Strojniški vestnik. Vse pravice pridržane. Prejeto v recenzijo: 2019-03-06 Prejeto popravljeno: 2019-04-03 Odobreno za objavo: 2019-04-11 Nelinearni model laserskega tvorjenja kapljic iz kovinske žice na osnovi dinamičnega ravnovesja dušenega vzmetnega nihala Andrej Jeromen* - Edvard Govekar Univerza v Ljubljani, Fakulteta za strojništvo, Slovenija V prispevku je predstavljen in ovrednoten nizko dimenzionalen matematični model procesa kontinuiranega laserskega tvorjenja kapljic iz kovinske žice, ki opiše bistveno dinamiko procesa, katere poznavanje je potrebno pri napovedovanju lastnosti kapljic in prediktivnem krmiljenju procesa. Pri procesu kontinuiranega laserskega tvorjenja kapljic podajamo kovinsko žico s konstantno hitrostjo v gorišče pulzirajočega anularnega laserskega žarka. Pri tem se konec žice pretali in oblikuje v rastočo visečo kapljico, ki se kasneje loči od žice. Opisani proces je zelo kompleksen zaradi svoje nelinearnosti ter nestacionarnosti in ker ga sestavlja preplet številnih fizikalnih pojavov. Eksperimentalne raziskave procesa so pokazale, da je v območju frekvence laserskih bliskov med 60 Hz in 190 Hz dinamika nastajajoče viseče kapljice zelo podobna dinamiki časovno odvisnega dušenega vzmetnega nihala ter da do ločitve kapljice od žice pride zaradi resonance med nihanjem viseče kapljice in frekvenco laserskih bliskov. Na osnovi teh eksperimentalnih rezultatov je bil formuliran matematični model viseče kapljice kot časovno spremenljive točkaste mase, viseče na vzmeti z dušilko. Vzbujanje zaradi periodičnega taljenja žice z laserskimi bliski je bilo opisano s spreminjanjem lege pritrdišča vzmeti. Model v obliki diferencialne enačbe 2. reda ima časovno odvisne koeficiente, ki so bili določeni na osnovi eksperimentalne analize rezultatov procesa. V modelu pride do ločitve kapljice od žice takrat, ko postane vsota vztrajnostne in gravitacijske sile večja od največje predvidene sile površinske napetosti, ki kapljico drži na koncu žice. Pri tem sem mora kapljica oddaljevati od žice, da se po ločitvi znova ne ujame na konec žice. Modelska enačba je bila rešena z numerično integracijo za vrednosti frekvence laserskih bliskov med 60 Hz in 190 Hz s korakom 10 Hz in za 31 enakomerno razporejenih vrednosti začetne mase viseče kapljice med 0,04*10-7 kg in 2,40*10-7 kg. Model je bil ovrednoten s primerjavo simuliranih in eksperimentalnih časovnih vrst vertikalne lege težišča viseče kapljice ter pripadajočih premerov ločenih kapljic pri različnih vrednostih frekvence laserskih bliskov. Primerjava je pokazala, da model dobro opiše bistveno dinamiko lege težišča viseče kapljice. Eksperimentalno je bilo potrjeno, da so zaznane razlike med simulirano in eksperimentalno določeno časovno vrsto posledica lateralnega nihanja viseče kapljice, ki ga model ne vključuje, in se vzbudi predvidoma zaradi odstopanj sistema od osne simetrije. Z modelom uspešno opišemo tudi rast in čas ločitve viseče kapljice kot posledico resonance nihanja viseče kapljice tako z osnovno frekvenco laserskih bliskov kot tudi z dvojno frekvenco laserskih bliskov. Posledično je bilo potrjeno dobro ujemanje odvisnosti modelskega in eksperimentalnega premera ločene kapljice od frekvence laserskih bliskov, vključno z bifurkacijo pri 120 Hz in nastalima dvema vejama vrednosti premera ločene kapljice pod 120 Hz. Odstopanje je opazno na mejah obravnavanega intervala frekvenc laserskih bliskov, kjer na ločevanje viseče kapljice poleg vertikalnih nihanj vpliva tudi lateralno nihanje kapljice in pojav Rayleigh-Plateaujeve nestabilnosti, ki ju model ne vključuje. V prispevku opisani model predstavlja eksperimentalno primerljiv nizko dimenzionalen matematični opis kompleksnega pojava kontinuiranega laserskega tvorjenja kapljic iz kovinske žice. Kot tak je zelo zanimiv za potencialno uporabo pri napovedovanju lastnosti tvorjenih kapljic in prediktivnem krmiljenju procesa. Vrednotenje rezultatov modela je med drugim pokazalo tudi na pomemben vpliv pojavov, kot sta Rayleigh-Plateaujeva nestabilnost in lateralno nihanje kapljice, ki lahko pri določenih frekvencah laserskih bliskov pomembno vplivata na izid procesa. Ključne besede: lasersko tvorjenje kapljic, sistem masa-vzmet-dušilka, nelinearni model, časovna vrsta, časovno-frekvenčna analiza *Naslov avtorja za dopisovanje: Univerza v Ljubljani, Fakulteta za strojništvo, Aškerčeva 6, 1000 Ljubljana, Slovenija, andrej.jeromen@fs.uni-lj.si SI 27 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, SI 28 © 2019 Strojniški vestnik. Vse pravice pridržane. Prejeto v recenzijo: 2019-02-20 Prejeto popravljeno: 2019-03-13 Odobreno za objavo: 2019-04-03 Eksperimentalna raziskava avtomobilske klimatske naprave s hladilnima sredstvoma R444a in R152a kot alternativo za R134a Mehmet Direk* - Mehmet Selguk Mert - Eren Soylu - Fikret Yuksel Univerza v Yalovi, Fakulteta za strojništvo, Turčija Uhajanje konvencionalnih hladilnih sredstev iz hladilnih sistemov v ozračje prispeva h globalnemu segrevanju. Hladilna sredstva so razvrščena glede na potencial globalnega segrevanja (GWP) in Evropska unija je v boju proti globalnemu segrevanju uvedla nekatere omejitve glede njihove uporabe in transporta. Kompresorski hladilni sistemi (VCR) so razširjeni na mnogih področjih, kamor spadajo tudi avtomobili. Večina avtomobilskih klimatskih naprav (AAC) uporablja hladilno sredstvo R134a z vrednostjo GWP 1300. Alternativna hladilna sredstva z manjšo vrednostjo GWP so npr. R1234ze(E), R152a in R444A. Sredstvo R444A vsebuje 83 % sredstva R1234ze(E), 12 % sredstva R32 in 5 % sredstva R152a, njegova vrednost GWP pa je 93 in tako izpolnjuje zahteve evropske zakonodaje. V predstavljeni študiji je bila preučena zmogljivost sistema AAC pri uporabi hladilnih sredstev R152a in R444A z nizko vrednostjo GWP kot alternativo za sredstvo R134a. Za določitev zmogljivosti in izboljšavo sistema je bil dodan tudi notranji prenosnik toplote (IHX). Uporabljen je bil eksperimentalni sistem AAC z ventili za odpiranje in zapiranje vodov za kapljevino in plin na IHX. V prvem delu so bila preizkušena hladilna sredstva R134a, R152a in R444A v sistemu AAC brez IHX, zasnovanem za sredstvo R134a. Nato je bil ovrednoten tudi vpliv IHX na delovanje sistema, napolnjenega s sredstvom R444A. Eksperimentalni sistem je bil zgrajen iz komponent avtomobilske klimatske naprave, vključno s ploščnim prenosnikom toplote med sesalnim vodom za kapljevino in plin. Sistem AAC je bil pred eksperimenti napolnjen s sredstvom R134a. Količina R134a je bila prilagojena tako, da je bilo doseženo najboljše hladilno število sistema. Za enako maso hladilnih sredstev je bila polnilna količina sredstev R152a in R444A določena glede na gostoto R134a v tekočem stanju. Količina polnjenja sredstev R134a, R152a in R444A je tako znašala 625 g, 460 g in 512 g. Deli sistema AAC so bili popolnoma izolirani pred zunanjimi vplivi. Nameščeni so bili merilni instrumenti, za prenos izmerjenih vrednosti pa je bil uporabljen sistem za zajemanje podatkov. Temperatura zračnih tokov na uparjalniku in kondenzatorju med obratovanjem sistema je znašala 27 °C oz. 35 °C. Vrtilna hitrost kompresorja se je povečevala od 750 vrt/min do 2750 vrt/min v korakih po 500 vrt/min. Meritve so bile opravljene po vzpostavitvi stacionarnega stanja. V prvem delu so bila preizkušena hladilna sredstva R134a, R152a in R444A brez prenosnika IHX. V drugem delu je bil nato aktiviran še IHX in opravljeni so bili eksperimenti s hladilnim sredstvom R444A. Rezultati so pokazali, da je bila največja hladilna zmogljivost dosežena s hladilnim sredstvom R152a, sledili pa sta sredstvi R134a in R444A. Hladilna zmogljivost sistema, napolnjenega s sredstvom R444A, je bila manjša kot pri R152a in R134a, ugotovljeno pa je bilo relativno povečanje zmogljivosti po aktiviranju prenosnika IHX. Največja vrednost hladilnega števila (COP) je bila za vse vrtilne hitrosti kompresorja ugotovljena pri sredstvu R152a, sledili pa sta mu sredstvi R134a in R444A. Hladilno število pri R444A se je povečalo po aktiviranju IHX. Najmanjša stopnja uničenja eksergije je bila pri vseh vrtilnih hitrostih kompresorja ugotovljena pri sredstvu R444A. Stopnja uničenja eksergije pri sredstvu R152a je bila nekoliko višja kot pri R444A, temu pa je sledil R134a. Ugotovljeno je bilo, da se več eksergije uniči pri aktiviranem prenosniku IHX in pri povišanih vrtilnih hitrostih kompresorja. Zaradi omejitev, ki jih postavlja zakonodaja EU, obstaja potreba po iskanju alternativ za hladilno sredstvo R134a. Rezultati študije bodo uporabni tudi za ekonomske analize. Ključne besede: R134a, R152a in R444A, mobilne klimatske naprave, COP SI 28 *Naslov avtorja za dopisovanje: Univerza v Yalovi, Fakulteta za strojništvo, Turčija, Oddelek za energetiko, Yalova, Turkey, mehmetdirek@hotmail.com Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, SI 29 © 2019 Strojniški vestnik. Vse pravice pridržane. Prejeto v recenzijo: 2018-08-29 Prejeto popravljeno: 2019-02-02 Odobreno za objavo: 2019-03-18 Metoda za ojačenje značilk v signalu vibracij pri vrtanju z večpasovnim spektralnim odštevanjem po valčni paketni dekompoziciji Youhang Zhou12* - Yong Li1 - Hanjiang Liu1 iUniverza v Xiangtanu, Šola za strojništvo, Kitajska 2Ministrstvo za šolstvo, Tehnično raziskovalno središče za kompleksne obdelovalne tehnologije in opremo, Kitajska Cilj študije je odprava težav pri izločanju značilk iz signalov vibracij, ki so obremenjeni z močnim šumom ozadja, konkretno odprava vpliva šuma pri analizi signala vibracij, zajetega med procesom vrtanja. Značilke v signalu vibracij pri vrtanju so bile ojačene z odpravo šuma po metodi valčne paketne dekompozicije (WPD) in spektralnega odštevanja. Predmet študije je odstranjevanje šuma iz signala in pridobivanje značilk za diagnosticiranje napak med obdelavo. Najprej so bili z dekompozicijo nadzornega signala in šuma v več plasti po metodi WPD pridobljeni signali v različnih frekvenčnih pasovih. Nato je bila izbrana velikost okvirja, primerna podfrekvenci signala šuma v različnih frekvenčnih pasovih, in pridobljena ocena povprečnega spektra signala v posameznih pasovih s kratko Fourierjevo transformacijo. Ustrezna podfrekvenca signala šuma je bila nato obdelana s spektralnim odštevanjem, podfrekvenca signala brez šuma pa je bila pridobljena z inverzno Fourierjevo transformacijo. Signali posameznih frekvenčnih pasov po spektralnem odštevanju so bili nato sešteti v čisti ojačeni signal. (1) Širokopasovni signal vibracij, obremenjen s šumom, je bil razdeljen v več podpasov po metodi WPD. S spektralnim odštevanjem in rekonstrukcijo signala je bil odpravljen vpliv šuma, ki je naložen na čisti signal. Na ta način je mogoče učinkovito ojačiti značilke v nadzornem signalu. (2) Z metodo WPMSS, ki zagotavlja osnovo za točnejšo identifikacijo in vrednotenje pri nadzoru procesa vrtanja, so bile ojačene značilke signala ob vstopu in izstopu svedra iz obdelovanca in opredeljena odvisnost med vrtalnim procesom in signalom vibracij. (3) Metoda WPMSS v primerjavi z originalnim spektralnim odštevanjem in valčno transformacijo učinkoviteje odstrani šum iz signala ter obenem ohrani šibke značilke vrtalnega procesa. Omejitev metode je v tem, da frekvenčna porazdelitev signalov, vključno s šibkimi značilkami in značilnimi frekvencami napak, zahteva delitev na frekvenčne pasove za točno oceno spektra. Prihodnje raziskave bodo osredotočene na kombiniranje različnih metod obdelave signalov za izboljšanje spektralnega odštevanja in na preučitev točnosti ocene spektra signala in šuma. Podan je predlog nove metode za ojačitev značilk signala vibracij pri vrtanju na osnovi večpasovnega spektralnega odštevanja po valčni paketni dekompoziciji. Metoda odpravlja znane težave pri izločanju šibkih značilk iz signala vibracij z močnim šumom ozadja. V primerjavi s tradicionalnimi metodami za odpravljanje šuma se je mogoče izogniti težavam, povezanim z izbiro pragov. Ključne besede: signal vibracij, zmanjšanje šuma v signalu, ojačenje značilk, valčna paketna dekompozicija, spektralno odštevanje *Naslov avtorja za dopisovanje: Univerza v Xiangtanu, Šola za strojništvo, Kitajska. zhouyouhang@xtu.edu.cn SI 29 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, SI 30 © 2019 Strojniški vestnik. Vse pravice pridržane. Prejeto v recenzijo: 2019-02-20 Prejeto popravljeno: 2019-03-13 Odobreno za objavo: 2019-04-03 Raziskava površinske hrapavosti delov iz jekel Hardox, obdelanih z abrazivnim vodnim curkom Alexandru Catalin Filip - Cristin Olimpiu Morariu - Laurentiu Aurel Mihail - Gheorghe Oancea* Transilvanska univerza v Brasovu, Romunija V zadnjih desetletjih beležimo visoko stopnjo rasti uporabe postopkov obdelave z abrazivnim vodnim curkom v industriji, predvsem zaradi njihove prilagodljivosti. Veljavnost rezultatov raziskav na tem področju pa je omejena na vrednosti parametrov procesa, materiale in specifične pogoje, ki so uporabljeni v eksperimentih. V članku je predstavljen znanstveni pristop k reševanju tega problema na področju obdelave jekel iz družine Hardox z abrazivnim vodnim curkom, saj predhodne raziskave rezultatov obdelave tega materiala niso bile zadostne. Raziskava je bila opravljena v okviru projekta industrijskega podjetja, ki izdeluje dele iz jekel Hardox po postopkih obdelave z abrazivnim vodnim curkom. Glavni cilj eksperimentalne raziskave je bil določitev odvisnosti med parametri procesa in kakovostjo površine delov iz jekel Hardox, obdelanih z abrazivnim vodnim curkom. V ta namen je bil razvit, analiziran in validiran matematični model odvisnosti. Izbrane so bile tri vhodne spremenljivke: vrsta materiala, podajalna hitrost rezanja in debelina dela. Določena je bila večfaktorska zasnova eksperimenta z dvema nivojema za vsakega od treh faktorjev. Vrednosti nivojev so bile določene na osnovi predhodnih raziskav in pragmatičnih potreb industrijskega podjetja. Vsak eksperiment je bil ponovljen trikrat za izboljšanje statistične stopnje zaupanja rezultatov. Po rezanju delov je bila analizirana površinska hrapavost na površini reza. Hrapavost je bila izmerjena s prenosnim merilnikom Mitutoyo SJ210 Surftest na sredini ravne ploskve, približno 4 mm od zgornjega in spodnjega roba. Za ovrednotenje površinske hrapavosti je bil uporabljen parameter Ra, ki je najbolj razširjen tako v industriji kot v raziskovalni sferi. Rezultati eksperimentov so bili ovrednoteni s statistično analizo srednjih vrednosti (ANOM) in s Pareto diagramom. Podajalna hitrost rezanja in debelina materiala pomembno vplivata na površinsko hrapavost, medtem ko je vrsta materiala statistično zanemarljiva. Podan je predlog matematičnih regresijskih odvisnosti za izračun parametra Ra za oba materiala, kakor tudi globalni model. Modeli so bili eksperimentalno validirani z zelo dobro stopnjo zaupanja. Regresijski modeli za napovedovanje vpliva parametrov so nelinearni. Eksperimentalna validacija je pokazala dobro stopnjo zaupanja globalnega modela za napovedovanje izhodnega parametra pri dveh materialih. Še višja stopnja zaupanja je dosežena s posebnim regresijskim modelom za vsak material posebej. Modeli so veljavni znotraj območja variabilnosti vhodnih parametrov, uporabljenih v eksperimentih. V prihodnjih raziskavah bodo razširjeni še z drugimi debelinami materiala, jekli Hardox in spremenljivimi vhodnimi parametri, kot je pretok abraziva. Praktična uporabnost novih matematičnih modelov je v tem, da bodo lahko operaterji strojev izbirali ustrezne vrednosti vhodnih parametrov za doseganje zahtevane kakovosti površin. Ključne besede: obdelava z abrazivnim vodnim curkom, jekla Hardox, površinska hrapavost, regresijski model SI 30 *Naslov avtorja za dopisovanje: Transilvanska univerza v Brasovu, Romunija, Oddelek za proizvodno strojništvo, Bdul Eroilor nr.29, Brasov, Romunija, gh.oancea@unitbv.ro Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, SI 31 © 2019 Strojniški vestnik. Vse pravice pridržane. Prejeto v recenzijo: 2018-12-03 Prejeto popravljeno: 2019-02-20 Odobreno za objavo: 2019-03-22 Ocena funkcionalnosti sistema za ustvarjanje mikroklime v stavbah v hladnih podnebjih s sončno energijo Karolis Januševičius - Juozas Bielskus - Vytautas Martinaitis Tehniška univerza Gediminas v Vilniusu, Fakulteta za okoljski inženiring, Litva Združeni narodi so kot enega od ciljev trajnostnega razvoja opredelili tudi oskrbo trgov s čisto in dostopno energijo. Toplotno ugodje je pomemben vidik učinkovite rabe energije, ki ima ključno vlogo pri zagotavljanju zdravja in dobrega počutja ljudi v stavbah. Največ energije v stavbah rabijo klimatski sistemi, katerih naloga je zagotavljanje toplotnega ugodja za uporabnike. Za zmanjšanje potrebne toplote za ogrevanje se uporabljajo strategije, kot je pasivno in aktivno izkoriščanje sončne energije ter toplotnih mas. Zagotavljanje toplotnega ugodja zahteva premišljene rešitve, ki slonijo na čistih, obnovljivih vire energije in ohranjajo druge vire. Za izpolnitev tega cilja so potrebni načini za povezovanje kriterijev toplotnega ugodja s projektnimi parametri, ki vplivajo na rabo energije in na okolje. Eden od pristopov k zmanjševanju vplivov zaradi rabe energije v trajnostnih zgradbah je z integracijo pasivnih in aktivnih tehnologij. Pokrivanje toplotnih potreb stavbe s solarnim ogrevanjem je omejeno z zmožnostmi akumulacije toplote. V članku je opisana metoda priprave diagramov za vrednotenje možnih zasnov v zgodnjih fazah projektiranja. Predlagani pristop lahko pomaga pri izbiri prave kombinacije projektnih parametrov na osnovi sezonskega toplotnega ugodja kot ciljnega parametra oz. glavne funkcije kompleksnih sistemov za ustvarjanje notranje mikroklime. Za preučevanje relacij med projektnimi spremenljivkami je bil uporabljen simulacijski model TRNSYS. Funkcionalnost sistema je bila ovrednotena z analizo ugodja v notranjih prostorih (po standardu EN ISO 7730). Simulacijski model upošteva toplotno konstanto notranjega prostora, sevalni tok skozi steklene površine, sončne kolektoije in aktivno akumulacijo, transport energije in izkoristek prenosa. Za poenostavitev procesa načrtovanja trajnostnih stavb je mogoče uporabiti diagrame za projektiranje. Ti diagrami so v predlaganem poteku dela pripravljeni s parametrično analizo v MATLAB-u in podajajo relacije med toplotnim ugodjem v ogrevalni sezoni ter projektnima parametroma deležem površine sončnih kolektorjev in toplotno maso. Funkcionalnost sistema je tukaj opredeljena s pričakovanim odstotkom nezadovoljnih (PPD), ki odraža sezonsko nezadovoljstvo uporabnikov na področju toplotnega ugodja. S tem pristopom je mogoče določiti projektne parametre, ki zagotavljajo sprejemljivo stopnjo zadovoljstva. Diagrami omogočajo opredelitev mejnih vrednosti za doseganje želenega toplotnega ugodja. Rezultati predstavljenega primera kažejo, da je mogoče pomanjkanje moči za ogrevanje prostorov (3/4 izračunane količine) kompenzirati s toplotno maso in s sončnimi kolektorji, ne da bi s tem ogrozili toplotno ugodje (vrednost PPD je pod 10 %). Delež površine sončnih kolektorjev (36 %) je zaradi večjih potreb največji pri najmanjši moči za ogrevanje prostorov (1/2 izračunane količine), toda ta možnost ne izpolnjuje kriterijev toplotnega ugodja (PPD > 17 %). Prednost predlagane metode je v tem, da omogoča hitro izbiro projektnih parametrov stavbe in sistema za ustvarjanje mikroklimatskih razmer že v zgodnji fazi projektiranja, omejitev metode pa je v računski zahtevnosti priprave diagramov za načrtovanje, ki se z vsako uvedbo novih spremenljivk še poveča. Diagrami prinašajo dodano vrednost v proces načrtovanja z normalizacijo projektnih spremenljivk in možnostmi ponovne uporabe. Algoritem za pripravo diagramov bi bilo mogoče še izboljšati z vključitvijo novih projektnih spremenljivk. Parametrično analizo bi bilo mogoče nadgraditi v iskanje kombinacij z večciljno optimizacijo in tako zmanjšati računsko zahtevnost priprave diagramov. Modeli toplotnega ugodja se razvijajo z vključevanjem metod na osnovi drugega zakona termodinamike in v prihodnje bo tako mogoče pripraviti tudi diagrame za zagotavljanje toplotnega ugodja z eksergijskimi kriteriji. Projektanti lahko s predlagano metodo snujejo kompleksne rešitve, ki upoštevajo več vidikov vpliva na raven toplotnega ugodja. Metoda vzpostavlja povezavo med različnimi projektnimi parametri (kot spremenljivkami) ter omogoča iskanje sinergij med različnimi aktivnimi in pasivnimi strategijami zagotavljanja toplotnega ugodja. Prednosti metode so: a) diagrame je treba pripraviti samo enkrat, medtem ko je znanje o medsebojnih odvisnostih trajno; b) metoda je primerna zamenjava za intuitivne metode »na palec« ter povečuje transparentnost v začetnih fazah načrtovanja; c) metoda omogoča iskanje kompromisov med raznimi pasivnimi in aktivnimi pristopi s sončnim ogrevanjem in akumulacijo toplote; d) trajnostno naravnani projektanti z njo dobijo orodje za izbiro in realizacijo trajnostnih kombinacij pasivne in aktivne rabe sončne energije z akumulacijo toplote. Ključne besede: sončno ogrevanje, orodje za zgodnje faze načrtovanja, toplotno ugodje *Naslov avtorja za dopisovanje: Tehniška univerza Gediminas v Vilniusu, Fakulteta za okoljski inženiring, Oddelek za energetiko v stavbah, Sauletekio av., 11, Vilnius, Litva, karolis.janusevicius @vgtu.lt SI 31 Strojniški vestnik - Journal of Mechanical Engineering 65(2019)4, SI 32 © 2019 Strojniški vestnik. Vse pravice pridržane. Prejeto v recenzijo: 2019-02-20 Prejeto popravljeno: 2019-03-13 Odobreno za objavo: 2019-04-03 Vpliv debeline prevleke na tok tekočine in trdnih delcev v centrifugalni črpalki Kaikai Luo12 - Yong Wang1* - Houlin Liu1 - Matevž Dular2 - Jie Chen1 - Zilong Zhang1 1 Univerza v Jiangsu, Raziskovalni center strojništva in tehnologije tekočin, Kitajska 2 Univerza v Ljubljani, Fakulteta za strojništvo, Slovenija Centrifugalne črpalke, ki delujejo v dvofaznem trdno-kapljevinskem toku, so nepogrešljive za transport snovi v metalurški, vodni in rudarski panogi. Zaradi posebnosti obratovalnih pogojev, mešanica trdne in tekoče faze povzroča poškodbe na delih stroja in s tem močno vpliva na stabilnost delovanja in na življenjsko dobo stroja. Navadno se za preprečitev tovrstnih poškodb uporabljajo premazni material z boljšo odpornostjo proti obrabi - na primer epoksidne smole, kompozitnimi najlonskimi premazi in poliuretanske prevleke. Po nanosu prevleke pa se spremenijo parametri stroja, kot sta debelina in kot lopatic, kar ima lahko velik vpliv na delovanje črpalke. V prispevku smo obravnavali vpliv premazov na tok v centrifugalni črpalki z nazivnim pretokom Qn = 20 m3/h, sesalno višino H = 22 m in izkoristkom n = 48 %. Rotor je bil prevlečen s poliuretanskim materialom debelin od 1 mm do 3 mm. V kombinaciji z numerično simulacijo smo podrobneje analizirali notranji pretok, nihanje tlaka in radialne sile na črpalko. Ujemanje napovedi integralnih karakteristik preko numerične simulacije z rezultati meritev je zelo dobro. Odstopanja med napovedjo vpliva in dejanskim izmerjenim zmanjšanjem tlačne višine in izkoristka so bila manjša od 5 %. Rezultati kažejo, da se s povečanjem debeline prevleke zmanjšata izkoristek in tlačna višina. Vpliv je sorazmeren z debelino prevleke. V primerjavi z originalnim (nezaščitenim rotorjem) se tlačna višina zmanjša za 3,74 % (1 mm debel nanos), 10,47 % (2 mm debel nanos) in za 24,24 % (3 mm debel nanos). Podobno se izkoristek zmanjša za 3,57 % (1 mm debel nanos), 6,22 % (2 mm debel nanos) oziroma 14,88 % (3 mm debel nanos). Izhodni tlak niha v nizkofrekvenčnem območju z jasnim vrhom, ki ustreza frekvenci prehoda lopatic. Ob obratovanju pri nazivnih pogojih se tlačne pulzacije z večanjem debeline nanosa povečujejo. Ob vstopu v rotor je prisotno večje območje z nizkim tlakom. Vstopni tlak se večanjem debeline nanosa najprej še poveča, nato pa zmanjša. Z naraščanjem debeline prevleke postane pretočnost v kanalu vse bolj nehomogena, posledično pa se z naraščanjem debeline prevleke povečajo tudi radialne sile na rotor. Ključne besede: debelina prevleke, dvofazna črpalka, numerična simulacija, nihanje tlaka, radialna sila SI 32 *Naslov avtorja za dopisovanje: Univerza v Jiangsu, Raziskovalni center strojništva in tehnologije tekočin, Zhenjiang, Kitajska, wylq@ujs.edu.cn Information for Authors All manuscripts must be in English. Pages should be numbered sequentially. The manuscript should be composed in accordance with the Article Template given above. The maximum length of contributions is 10 pages. Longer contributions will only be accepted if authors provide juastification in a cover letter. For full instructions see the Information for Authors section on the journal's website: http://en.sv-jme.eu . SUBMISSION: Submission to SV-JME is made with the implicit understanding that neither the manuscript nor the essence of its content has been published previously either in whole or in part and that it is not being considered for publication elsewhere. All the listed authors should have agreed on the content and the corresponding (submitting) author is responsible for having ensured that this agreement has been reached. The acceptance of an article is based entirely on its scientific merit, as judged by peer review. Scientific articles comprising simulations only will not be accepted for publication; simulations must be accompanied by experimental results carried out to confirm or deny the accuracy of the simulation. Every manuscript submitted to the SV-JME undergoes a peer-review process. The authors are kindly invited to submit the paper through our web site: http://ojs.sv-jme.eu. The Author is able to track the submission through the editorial process - as well as participate in the copyediting and proofreading of submissions accepted for publication - by logging in, and using the username and password provided. SUBMISSION CONTENT: The typical submission material consists of: - A manuscript (A PDF file, with title, all authors with affiliations, abstract, keywords, highlights, inserted figures and tables and references), - Supplementary files: • a manuscript in a WORD file format • a cover letter (please see instructions for composing the cover letter) • a ZIP file containing figures in high resolution in one of the graphical formats (please see instructions for preparing the figure files) • possible appendicies (optional), cover materials, video materials, etc. Incomplete or improperly prepared submissions will be rejected with explanatory comments provided. In this case we will kindly ask the authors to carefully read the Information for Authors and to resubmit their manuscripts taking into consideration our comments. COVER LETTER INSTRUCTIONS: Please add a cover letter stating the following information about the submitted paper: 1. Paper title, list of authors and their affiliations. 2. Type of paper: original scientific paper (1.01), review scientific paper (1.02) or short scientific paper (1.03). 3. A declaration that neither the manuscript nor the essence of its content has been published in whole or in part previously and that it is not being considered for publication elsewhere. 4. State the value of the paper or its practical, theoretical and scientific implications. What is new in the paper with respect to the state-of-the-art in the published papers? Do not repeat the content of your abstract for this purpose. 5. We kindly ask you to suggest at least two reviewers for your paper and give us their names, their full affiliation and contact information, and their scientific research interest. The suggested reviewers should have at least two relevant references (with an impact factor) to the scientific field concerned; they should not be from the same country as the authors and should have no close connection with the authors. FORMAT OF THE MANUSCRIPT: The manuscript should be composed in accordance with the Article Template. The manuscript should be written in the following format: - A Title that adequately describes the content of the manuscript. - A list of Authors and their affiliations. - An Abstract that should not exceed 250 words. The Abstract should state the principal objectives and the scope of the investigation, as well as the methodology employed. It should summarize the results and state the principal conclusions. - 4 to 6 significant key words should follow the abstract to aid indexing. - 4 to 6 highlights; a short collection of bullet points that convey the core findings and provide readers with a quick textual overview of the article. These four to six bullet points should describe the essence of the research (e.g. results or conclusions) and highlight what is distinctive about it. - An Introduction that should provide a review of recent literature and sufficient background information to allow the results of the article to be understood and evaluated. - A Methods section detailing the theoretical or experimental methods used. - An Experimental section that should provide details of the experimental set-up and the methods used to obtain the results. - A Results section that should clearly and concisely present the data, using figures and tables where appropriate. - A Discussion section that should describe the relationships and generalizations shown by the results and discuss the significance of the results, making comparisons with previously published work. (It may be appropriate to combine the Results and Discussion sections into a single section to improve clarity.) - A Conclusions section that should present one or more conclusions drawn from the results and subsequent discussion and should not duplicate the Abstract. - Acknowledgement (optional) of collaboration or preparation assistance may be included. Please note the source of funding for the research. - Nomenclature (optional). Papers with many symbols should have a nomenclature that defines all symbols with units, inserted above the references. If one is used, it must contain all the symbols used in the manuscript and the definitions should not be repeated in the text. In all cases, identify the symbols used if they are not widely recognized in the profession. Define acronyms in the text, not in the nomenclature. - References must be cited consecutively in the text using square brackets [1] and collected together in a reference list at the end of the manuscript. - Appendix(-icies) if any. SPECIAL NOTES Units: The SI system of units for nomenclature, symbols and abbreviations should be followed closely. Symbols for physical quantities in the text should be written in italics (e.g. v, T, n, etc.). Symbols for units that consist of letters should be in plain text (e.g. ms-1, K, min, mm, etc.). Please also see: http://physics.nist.gov/cuu/pdf/sp811.pdf . Abbreviations should be spelt out in full on first appearance followed by the abbreviation in parentheses, e.g. variable time geometry (VTG). The meaning of symbols and units belonging to symbols should be explained in each case or cited in a nomenclature section at the end of the manuscript before the References. Figures (figures, graphs, illustrations digital images, photographs) must be cited in consecutive numerical order in the text and referred to in both the text and the captions as Fig. 1, Fig. 2, etc. Figures should be prepared without borders and on white grounding and should be sent separately in their original formats. If a figure is composed of several parts, please mark each part with a), b), c), etc. and provide an explanation for each part in Figure caption. The caption should be self-explanatory. Letters and numbers should be readable (Arial or Times New Roman, min 6 pt with equal sizes and fonts in all figures). Graphics (submitted as supplementary files) may be exported in resolution good enough for printing (min. 300 dpi) in any common format, e.g. TIFF, BMP or JPG, PDF and should be named Fig1.jpg, Fig2.tif, etc. However, graphs and line drawings should be prepared as vector images, e.g. CDR, AI. Multi-curve graphs should have individual curves marked with a symbol or otherwise provide distinguishing differences using, for example, different thicknesses or dashing. Tables should carry separate titles and must be numbered in consecutive numerical order in the text and referred to in both the text and the captions as Table 1, Table 2, etc. In addition to the physical quantities, such as t (in italics), the units [s] (normal text) should be added in square brackets. Tables should not duplicate data found elsewhere in the manuscript. Tables should be prepared using a table editor and not inserted as a graphic. REFERENCES: A reference list must be included using the following information as a guide. Only cited text references are to be included. Each reference is to be referred to in the text by a number enclosed in a square bracket (i.e. [3] or [2] to [4] for more references; do not combine more than 3 references, explain each). No reference to the author is necessary. References must be numbered and ordered according to where they are first mentioned in the paper, not alphabetically. All references must be complete and accurate. Please add DOI code when available. Examples follow. Journal Papers: Surname 1, Initials, Surname 2, Initials (year). Title. Journal, volume, number, pages, DOI code. [1] Hackenschmidt, R., Alber-Laukant, B., Rieg, F. (2010). Simulating nonlinear materials under centrifugal forces by using intelligent cross-linked simulations. Strojniški vestnih - Journal of Mechanical Engineering, vol. 57, no. 7-8, p. 531-538, DOI:10.5545/sv-jme.2011.013. Journal titles should not be abbreviated. Note that journal title is set in italics. Books: Surname 1, Initials, Surname 2, Initials (year). Title. Publisher, place of publication. [2] Groover, M.P. (2007). Fundamentals of Modern Manufacturing. John Wiley & Sons, Hoboken. Note that the title of the book is italicized. Chapters in Books: Surname 1, Initials, Surname 2, Initials (year). Chapter title. Editor(s) of book, book title. Publisher, place of publication, pages. [3] Carbone, G., Ceccarelli, M. (2005). Legged robotic systems. Kordic, V., Lazinica, A., Merdan, M. (Eds.), Cutting Edge Robotics. Pro literatur Verlag, Mammendorf, p. 553576. Proceedings Papers: Surname 1, Initials, Surname 2, Initials (year). Paper title. Proceedings title, pages. [4] Štefanic, N., Martinčevic-Mikic, S., Tošanovic, N. (2009). Applied lean system in process industry. MOTSP Conference Proceedings, p. 422-427. Standards: Standard-Code (year). Title. Organisation. Place. [5] ISO/DIS 16000-6.2:2002. Indoor Air — Part 6: Determination of Volatile Organic Compounds in Indoor and Chamber Air by Active Sampling on TENAX TA Sorbent, Thermal Desorption and Gas Chromatography using MSD/FID. International Organization for Standardization. Geneva. WWW pages: Surname, Initials or Company name. Title, from http://address, date of access. [6] Rockwell Automation. Arena, from http://www.arenastmuiation.com, accessed on 200909-07. EXTENDED ABSTRACT: When the paper is accepted for publishing, the authors will be requested to send an extended abstract (approx. one A4 page or 3500 to 4000 characters). The instruction for composing the extended abstract are published on-line: http://www.sv-ime.eu/mformation-for-authors/ . COPYRIGHT: Authors submitting a manuscript do so on the understanding that the work has not been published before, is not being considered for publication elsewhere and has been read and approved by all authors. The submission of the manuscript by the authors means that the authors automatically agree to transfer copyright to SV-JME when the manuscript is accepted for publication. All accepted manuscripts must be accompanied by a Copyright Transfer Agreement, which should be sent to the editor. The work should be original work by the authors and not be published elsewhere in any language without the written consent of the publisher. The proof will be sent to the author showing the final layout of the article. Proof correction must be minimal and executed quickly. Thus it is essential that manuscripts are accurate when submitted. Authors can track the status of their accepted articles on http://en.sv-ime.eu/. PUBLICATION FEE: Authors will be asked to pay a publication fee for each article prior to the article appearing in the journal. However, this fee only needs to be paid after the article has been accepted for publishing. The fee is 380 EUR (for articles with maximum of 6 pages), 470 EUR (for articles with maximum of 10 pages), plus 50 EUR for each additional page. The additional cost for a color page is 90.00 EUR. These fees do not include tax. Strojniški vestnik -Journal of Mechanical Engineering Aškerčeva 6, 1000 Ljubljana, Slovenia, e-mail: info@sv-jme.eu http://www.sv-jme.eu Contents Papers Andrej Jeromen, Edvard Govekar: Nonlinear Dynamic Force Balance Mass-spring-damper Model of Laser Droplet Generation from a Metal Wire Mehmet Direk, Mehmet Selcuk Mert, Eren Soylu, Fikret Yüksel: Experimental Investigation of an Automotive Air Conditioning System Using R444A and R152a Refrigerants as Alternatives of R134a Youhang Zhou, Yong Li, Hanjiang Liu: Feature Enhancement Method for Drilling Vibration Signals by Using Wavelet Packet Multi-band Spectral Subtraction Alexandru Catalin Filip, Cristin Olimpiu Morariu, Laurentiu Aurel Mihail Gheorghe Oancea: Research on the Surface Roughness of Hardox Steel Parts Machined with an Abrasive Waterjet Karolis Januševičius, Juozas Bielskus, Vytautas Martinaitis: Functionality Assessment of Building a Microclimate System Utilising Solar Energy in a Cold Climate Kaikai Luo, Yong Wang, Houlin Liu, Matevž Dular, Jie Chen, Zilong Zhang Effect of Coating Thickness on a Solid-Liquid Two-Phase Flow Centrifugal Pump under Water Medium 9770039248001