Swing check valve

IN

DEGREE PROJECT MECHANICAL ENGINEERING, SECOND CYCLE, 30 CREDITS

STOCKHOLM SWEDEN 2017 ,

Swing check valve

characterization: 3D CFD

validation of one dimensional models used in RELAP5

GIACOMO VAVASSORI

Abstract

In a previous thesis work a swing check valve was studied with CFD anal- ysis in order to find correlations that could provide a good input for a one- dimensional model of the same. In this document, starting from the previous thesis results and using the model by Li and Liou as the reference work, a check valve is investigated and the hydraulic torque coefficients identified. In this way it becomes possible to analyze the behavior of the same valve with a 1D code called RELAP5.

The first part of the work was dedicated to understanding the dynamics lying behind the movement of swing check vale, and to the construction of a suitable 3D CFD model being able to find the required coefficients. The results were subsequently elaborated and implemented in the RELAP5 code, in order to run the 1D simulations. In the second part of the job, several transient simulations of different pipelines were conducted with the 1D model, monitoring in particular the closure of the valve over time.

In the end, the data obtained in RELAP5 were compared to those from equivalent 3D CFD analysis and from an alternative 1D approach by Adamkowski.

Although not completely matching, the results showed that the 1D model by

Li and Liou has a good ability to accurately simulate the valve. Further work

could be surely done on particular topics in order to improve and tune this

model, making it a consistent and cheaper alternative to 3D simulations, able

to accurately simulate a wide range of cases.

Acknowledgement

I would like to thank Ori Levin for assisting me especially in the last and most important stages of this work, and for allowing me to conduct this thesis at FS Dynamics Sweden. I would also like to say thanks to Elin and Emil who helped me throughout the many phases of my work, especially when things were not going as supposed.

I would like to thank my supervisor and examiner at KTH Luca Brandt that gave me useful tips and helped me during the thesis writing.

Last but not least I want to say thanks to my family that supported me

throughout all my studies and made it possible for me to finally become an

engineer.

Contents

1 Introduction 7

1.1 Thesis Goal . . . . 7

1.2 Swing check valve characterization . . . . 8

1.2.1 Previous Works . . . . 9

2 Background 11 2.1 Computational Fluid Dynamics . . . . 11

2.1.1 Turbulence Models . . . . 13

2.2 Model for the Swing Check Valve Dynamics . . . . 14

2.2.1 Calculating Hydraulic Torque Coefficients According to Li and Liou’s Model . . . . 17

2.2.2 Alternative 1D Model for Hydraulic Torque . . . . 19

2.3 About RELAP5 . . . . 20

3 Numerical Set-Up 22 3.1 Geometry . . . . 22

3.2 Boundary Conditions . . . . 24

3.3 Convergence Criteria . . . . 25

3.4 Sensitivity Analysis . . . . 26

3.5 STAR-CCM+ Model Set-Up . . . . 26

3.5.1 Mesh Generation . . . . 26

3.5.2 Fluid - Body Interaction . . . . 29

3.6 RELAP5 Model Set-Up . . . . 29

3.6.1 Check Valve Modeling Input . . . . 29

3.6.2 System Description . . . . 31

4 Results 33 4.1 Valve Investigation with 3D CFD . . . . 33

4.1.1 Sensitivity Analysis For Steady State Simulations . . . . . 33

4.1.2 Steady State Analysis: Data Elaboration and Results . . 41

4.1.3 Sensitivity Analysis For Transient Simulations . . . . 45

4.1.4 Transient Simulations: Data Elaboration and Results . . 47

4.1.5 Valve Geometry Scaling . . . . 52

4.2 RELAP5 Simulations . . . . 54

4.3 RELAP5 Models versus STAR-CCM+ . . . . 56

5 Summary and Conclusions 60

5.1 Future Works . . . . 62

Appendices 65

A Matlab Script for C _HR Regressions 66

B RELAP5 Input for Sr Model 69

C RELAP5 Valve Definition for Cm-Cq Model 92

List of Tables

3.1 Fluid Operative Conditions . . . . 23

3.2 Valve Door Characteristics . . . . 23

3.3 Mesh for Sensitivity and Steady State Analysis . . . . 27

3.4 Final Mesh for Unsteady State Analysis . . . . 29

4.1 Mesh Sensitivity Analysis . . . . 35

4.2 Steady State Sensitivity Analysis . . . . 36

4.3 Solver Sensitivity Analysis . . . . 36

4.4 Roughness effects on the system . . . . 37

4.5 Turbulence models comparison . . . . 38

4.6 C _HS variation with different Re . . . . 38

4.7 Mesh Tested . . . . 45

4.8 Boundary Conditions for Simulated Cases . . . . 47

4.9 Flow Accelerations for RELAP5 Runs . . . . 55

4.10 Differences in Normalized Angle at Norm. Time = 0.6 . . . . 57

4.11 Door Average Velocities and Accelerations for Different Models . 59

List of Figures

2.1 Schematic Representation of the Swing Check Valve . . . . 16

3.1 Geometry of the Valve . . . . 23

3.2 Angle Reference System . . . . 31

3.3 Scheme for the First Pipeline . . . . 32

3.4 Scheme for the Second Pipeline . . . . 32

4.1 Sectionf of the valve at 35 ^◦ opening . . . . 34

4.2 Normal Probabilty Plot for Forward CHS . . . . 42

4.3 Collected Data and Trendlines for Static Coefficient . . . . 43

4.4 Computed Values for C _M and C _Q Coefficients . . . . 44

4.5 Inlet Flow Velocity for Case 4 Simulation . . . . 48

4.6 Interpolated C _HR at 10 rad/s . . . . 49

4.7 Sample from regression result for C _HR . . . . 50

4.8 R ² value for each constant-angle regression . . . . 50

4.9 Torques for scaled and normal-sized valves . . . . 53

4.10 Rotational Coefficient for scaled and normal-sized valves . . . . . 53

4.11 Example of Mass Flow Rate Over Time for -140 m/s ² case . . . 54

4.12 Valve Door Closure over Time for the Different Runs . . . . 55

4.13 Comparison between 1D and 3D Models . . . . 56

4.14 Comparison between 1D and 3D Models with Normalized Times 58 4.15 Valve Closure for Pump Trip Comparison . . . . 59

5.1 Time Required for Valve Closure between 1D and 3D Model . . . 61

Abbreviation Description

CFD Computational Fluid Dynamics

S R Li and Liou’s model as referred in RELAP5 C M − C Q Adamkowski’s model as referred in RELAP5

θ Valve opening angle [deg]

ω Angular velocity [rad/s]

θ ˙ Angular velocity [rad/s]

˙

ω Angular Acceleration [rad/s ² ]

θ ¨ Angular Acceleration [rad/s ² ]

ρ Density [kg/m ³ ]

A v Valve disc frontal area [m ² ]

Re Reynolds number

µ Dynamic Viscosity [P a · s]

L Lever arm length [m]

c Speed of sound [m/s]

C HS Stationary hydraulic torque coefficient

C _HR Rotational hydraulic torque coefficient

Chapter 1

Introduction

1.1 Thesis Goal

An incompressible fluid, flowing in a rigid pipe, can eventually generate a so- called water hammer. This phenomenon occurs when a flow is rapidly stopped, the inertia of the fluid column itself will cause a pressure surge that could even result in catastrophic effects. Check valves are components meant to prevent the flow from streaming back along a pipe. Some possible reasons for a fluid to start flowing backwards are for example presence of any kind of pipe break or a pump shut down transient. Check valves are therefore meant to immediately stop the flow before the reverse flow happen, this in order to avoid the water hammer effects [1][2].

There are some common valve system designs able to interrupt the flow when

a pressure drop occurs. These components are necessary to ensure stability

in the overall system when for example two or more pumps are operating in

parallel, and one needs to be shut off. Several valve types exist; the most

recurrent are for example the swing check valve, the tilting disc, the double disc

and the lift valve. All these operate in a similar way, thanks to their design

which is meant to exploit the force exerted by the fluid onto the valve doors in

order to maintain them opened in the desired flow conditions, while at the same

time being able to use the same fluid force to shut off the valve with adverse

stream flows. In this master thesis work the focus is posed on swing check valve

mechanisms, which are particularly prone to generate water hammers due to the

fact that their simple design relies on the force generated by the reverse flow for

getting the valve door completely closed. The valve is indeed shut down quite

rapidly once the backwards flow has initiated, but its subsequent almost instant

disruption creates the perfect conditions for a water hammer. Water slam can

however be prevented if the valve closes at the same time as when the flow

reaches a 0 m/s velocity. Different design of check valves are capable of rapidly

interrupt flow before slam occurs thanks to springs or lever externally connected

to the valve door. In [3], several tests were run on different valve designs in

order to measure the back-flow velocity as a function of the flow deceleration.

Swing check valves were placed among the worst performing devices in terms of preventing slam under high flow deceleration.

In order to be able to estimate in advance the water hammer intensity, it is of extreme importance to be able to fully characterize the behavior of the check valve during a flow transient. To prevent water hammers, a check valve should be able to either rapidly close before appreciable back-flow occurs, or close very slowly afterward [2].

Swing check valves design consists of a disc rotating around a fixed axis. It is therefore sufficient to solve the Newton’s second law for momentum to fully describe the disc movement:

T tot = I tot θ ¨ (1.1)

Many studies have been performed on swing check valves in order to find a systematic way to fully determine their behavior under normal operative con- ditions in a pipe circuit system. However, while a static modelling has been successful thanks to the fact that it is relatively easy to determine flow condi- tions in a steady state flow, the behavior of such components under dynamics conditions is still an interesting topic and is open for further investigations.

When designing a pipeline it is indeed hard to know in advance what will be the possible negative acceleration that the system is going to witness. Using a computing software will surely facilitate the designing process, but the solution provided will not always be accurate or easy to get. One software that can be used is RELAP5, which is a one dimensional code that determines pipe systems’

interesting values. However this code was performing poorly when modeling transient flows across swing check valves. In order to improve the program reliability different alternative hydraulic torque models were implemented in the code, and the goal of this thesis work is indeed to test the validity of one of these approaches through a comparison of results obtained between 3D CFD simulation with STAR-CCM+ and 1D computations with RELAP5.

1.2 Swing check valve characterization

As mentioned before, pressure surges have a negative impact on pipes, they generate stresses that can result in ruptures and major damages if the pressure is high enough. The first engineers to thoroughly investigate this phenomenon were Joukowski and Allievi, which came to similar results independently. Their equation points out that pressure surges are a function of the flow velocity; the higher the velocity, the greater the pressure wave generated. For a swing check valve the velocity of interest is the backwards flow fluid speed, which should ideally be equal to zero in order to avoid any surge, but this happens seldom.

The common formulation for the Joukowski equation is

∆P = ρc∆V (1.2)

The main obstacle to a complete valve modeling is the fact that it is hard to accurately determine the various torque components acting on the valve door (from eq. 1.1), which are needed to solve the governing equation for swinging disc. These components can generally be divided into the following [5]:

T tot = T W + T f + T e + T h (1.3) where T _W is the torque due to the submerged weight of the disc, T _f is the frictional torque component acting on the rotational shaft, T _e is an eventual externally applied torque and T _h is the moment due to fluid flowing, hereby identified as the hydraulic torque. And it is the hydraulic torque component in particular to represent the most difficult element to characterize for a successful and easily accessible physical model. It should be determined by knowing the pressure on every point of the disc surface. But, since it is clearly impossible to directly measure the infinitesimal pressure in any kind of laboratory experiment, the easiest alternative approach commonly used nowadays to characterize the valve consists of measuring the pressure difference and the flow rate across the door. This leads to a so-called flow coefficient C v [4] which is part of the original numerical model for check valve characterization used in the aforementioned RELAP5 code. But this approach simply based on pressure difference over the valve is not enough accurate.

In the past years many alternative models were proposed in order to find a solution to this problem. Of particular interest is the method provided by Li and Liou for transient characteristic of check valves [5], based also on a previous work by Ellis and Mualla [6]. The main idea was that the hydraulic component of the torque acting on the disc could be divided in two parts, one acting on a fixed disc and subsequently called stationary component T _HS , and one accounting for the changing pressure on the swinging disc, hereby defined as the rotational component T HR . The two components could be calculated through means of two coefficients, respectively C HS and C HR , that were experimentally determined, and varying for each valve design. This model was tested experimentally and proved valid and accurate for relatively small flow accelerations, being capable of reliably predict closing time for the tested check valve. This new model has also been implemented in RELAP5 sotware, and it is now available as an alternative approach for 1D check valve characterization. However, the program requires the two hydraulic torque coefficients as an input, and almost no valve manufacturer provides them.

Using a computational fluid dynamic software makes it possible to estimate these value for any swing check valve. The advantage of performing a compu- tational experiment relies on the fact that the hydraulic torque can be easily computed from the effective pressure acting on the valve door surface, while it is clearly impossible to do that for a physical experiment in a laboratory.

1.2.1 Previous Works

There are already some extensive works which have been focusing on CFD in-

vestigation for check valves during transient flows. In particular Turesson [7]

performed a comparison between three one-dimensional models adopted in RE- LAP5 and the 3D CFD simulation for the same valve, pointing out the already known issues concerning the mono-dimensional modeling methodologies. Jans- son and L¨ ovmark [9] also studied a check valve with the final goal of comparing the computed data with physical experimentation, however the work is limited only to the CFD analysis.

Boqvist [8] studied a valve with CFD simulation, aiming to find trends and correlation in dynamic valve behavior, which could be used to improve or im- plement the models adopted in RELAP5. Eriksson [19] presented an analysis of the same down-scaled check valve under a pump trip event.

This thesis work will be based on the aforementioned documents, which are

also used as a comparison and a guideline for some parts of the present master

thesis. The investigated valve geometry is indeed the same as in [8] [19], and

the final goal of the present work is partly based on suggested future works

contained in those papers.

Chapter 2

Background

In this chapter the theory behind check valve dynamics and computational fluid dynamics will be explained, as well as a description of the principal issues to consider when setting up the simulation. In addition to that, a brief illustration of the RELAP5 code will be given.

According to [5], a better reliability compared to the C v method can be achieved with the new model, and this is why it has been added into RELAP5.

With the CFD code STAR-CCM+ a simulation of the swing check valve will be performed under different flow conditions in order to firstly find out the needed coefficients. These will be used afterwards in the one dimensional software to compute the pipe system interesing parameters for different cases involving check valve motion, and then, finally, the same cases will be simulated with the CFD software, in order to check the pros and cons of the two different methods.

2.1 Computational Fluid Dynamics

Computational Fluid Dynamics is based on the fundamental conservation equa- tions for mass, momentum and energy::

∂ρ

∂t + ∇ · ρV
= 0 (2.1)

The momentum conservation equation is (Cauchy):

∂

∂t (ρV) + ∇ · (ρVV) = ρg + ∇σ ij (2.2) And the energy conservation is expressed as:

∂(ρi)

∂t + ∇ · (ρiu) = −p∇u + ∇ · (k∇T ) + Φ (2.3)

where i is the internal thermal energy of the fluid, k is the thermal conductivity

and u is the total energy as a sum of thermal, kinetic and potential energy.

In order to solve the momentum equation for a velocity field, additional equations called constitutive equations have to be provided. There are several of these which characterize physical properties such as fluid viscosity, density etc., the most important is the equation of state which has several expression depending on the fluid modeled:

p = p(ρ, T ) i = i(ρ, T ) (2.4)

For Newtonian incompressible fluids the Cauchy equation 2.2 can be used to formulate the so-called Navier-Stokes Equations, which are a fundamental result for fluid mechanics. For an incompressible flow case, they can be expressed in a single vector equation [14]:

ρ DV

Dt = −∇P + ρg + µ∇ ² V (2.5)

The Navier-Stokes equations allow to completely determine the fluid behav- ior but cannot be analytically solved, if not for just very specific and simple cases. Alternative solutions to these equations have been proposed and investi- gated.

The Direct Numerical Solution of the equations is theoretically possible, but in reality it is not convenient at all. DNS involves attempting to solve the unsteady flow for every scale level for both space and time [14]. This is usually a quite wide range, especially in normal industrial application such as the case investigated in the present work. Indeed, the difference in scale of the flow phenomena becomes higher with increasing Reynolds number. This involves that a very small time step and a very 3D refined grid must be used to compute the solution in the whole domain, with a substantial impact on the computational power required. With the calculation power available today it is not possible at all to recur to this way to solve common problem.

Since in many case it is not necessary to achieve a very detailed solution for the flow, a more convenient way to model turbulent flows has been developed;

the solution consists of using an alternative formulation of the Navier-Stokes equations, based on the possibility to express velocities and pressure as a sum of a time-averaged component and a fluctuating one, this procedure is called Reynolds decomposition and it is generically expressed by:

φ(t) = ¯ φ + φ ⁰ (t) (2.6)

The so-called RANS (Reynolds Averaged Navier-Stokes) equations can be therefore derived from the previous equations.

The formulation of RANS equations introduces at the same time additional requirements, additional equations are indeed necessary to close the system.

These can be found in the various turbulence models proposed and discussed in

the next section.

2.1.1 Turbulence Models

Many of the nowadays used turbulence models are based on the assumption that an analogy exists between the Reynold stresses acting in the mean flow and the viscous stresses, in this way it is possible to model the Reynolds stresses directly as a function of mean flow properties. In addition to that, a distinction among models can be done according to the number of transport equations introduced.

There are zero, one, two and seven equations model.

Zero Equation Model

The Mixing Length model is the only zero model introduced here and it is based on the idea that the dynamic turbulent viscosity can be expressed as a product of a turbulent velocity scale θ and a turbulent length scale l:

µ _t = Cρθl (2.7)

This model is useful for 2D flows with smooth changes in flow direction (no separation) and where the turbulence creation and dissipation are balanced. If this condition is true, turbulence properties are proportional to the mean flow properties (the mixing length scale l is proportional to the mean flow length scale L) and they can be described with simple formulas.

One Equation Model

The Spalart-Allmaras model belongs to the One Model class, similarly to the Mixing Length, it can be used mostly when there is almost no separation and indeed it was derived mainly for aerospace applications. This model works fine with low Reynolds number but does not provide good results for typical industrial application cases, like the one investigated in the present work. The eddy viscosity is calculated through equation:

µ t = ρf v ν ˜ (2.8)

The transport equation introduced by this model is [12]:

∂

∂t ρ˜ ν
+ ∇ · (ρ˜ ν ¯ U) = 1

σ ∇ · [(µ + ρ˜ ν)∇˜ ν] + P ν ˜ + S ν ˜ (2.9) Where σ is a model coefficient, and P ν ˜ and S ν ˜ are the source and production terms.

Two Equations Model

This category contains the most used model for common applications where the

flow presents high Reynolds number and separation can occur. The first and

well-known model is the k − ε which presents good compromises between ro-

bustness and accuracy. This model has been used widely and several forms have

been developed to best suit the simulated case properties in terms of Reynolds

number and wall treatment. The k − ε introduces two additional transport equation for turbulent kinetic energy k and turbulent dissipation rate ε [13]:

∂ρk

∂t + ∇ · (ρkU) = ∇ ·  µ t

σ k

∇k



+ 2µ t S ij · S ij − ρε (2.10)

∂ρε

∂t + ∇ · (ρεU) = ∇ ·  µ _t σ ε

∇ε

 + C 1ε

ε

k 2µ t S ij · S ij − C 2ε ρ ε ²

k (2.11)

These two equations contains 5 coefficients which vary for the different mod- els adopted. The k − ε is the most adopted solution and this is due to its ability to provide reliable results with good performances, especially for internal flow problems.

Alternatively to this one, there is a second two-equations model called SST k−

ω, which is based on turbulent kinetic energy as well, but it takes into account the specific dissipation rate ω = ε/k instead of the turbulent dissipation rate.

The main advantage of this model compared to the k − ε relies on its improved performance while describing the flow in the boundary layers, especially with strong pressure gradients. Compared to the original version of the k − ω, it is also less sensitive to the initial arbitrary values for the free stream [13], this makes it a good alternative for general CFD applications with internal flow problems. The new two transport equations introduced by this model are:

∂ρk

∂t + ∇ · (ρkU) = ∇ ·



µ + µ _t σ k



∇k



+ P _k − β ^∗ ρkω (2.12)

∂ρω

∂t + ∇ · (ρωU) = ∇ ·



µ + µ t

σ _ω1



∇ω

 + γ 2



2ρS ij · S ij − 2 3 ρω ∂U i

∂x _j δ ij

 +

− β 2 ρω ² + 2 ρ σ _ω2 ω

∂k

∂x _k

∂ω

∂x _k

(2.13)

2.2 Model for the Swing Check Valve Dynamics

A swing check valve consists essentially of a disc pivoting around a constrained rotational axis, which is therefore the only degree of freedom allowed for the valve door. Before proceeding, it is important to set a common reference sys- tem which will be used throughout the whole report. The momentum will be considered positive if exerting a counterclockwise rotation on the valve door. As mentioned in section 1.1, the motion can be fully expressed by the second law of dynamics for rotation (eq. 1.1).

Even if this equation can totally describe the movement, it is not really enough detailed and does not give an idea of which are the causes that originates the motion. A better formulation from [5] is:

T W + T f + T e + T H = (I + I e + I f ) ˙ ω (2.14)

Before analyzing eq. 2.14, it is important to introduce the concept of the center of mass, which is a reference point where the weight force of the disc is applied. This point is of extreme importance as it is used for calculating the angle that affects the weight torque magnitude as well as the hydraulic one. It is also used to calculate the disc inertia, which plays an important role in the equilibrium equation, as explained later on.

A check valve operates constantly in contact with a fluid, either flowing or not, which exerts a force on the door structure. The effect due to the fluid will be denoted as the hydraulic torque, that is the subject of the investigation. Hy- draulic torque can be physically described as the result of pressure and viscous effects. The stresses acting on the surface can be expressed by the symmetric stress tensor, or Cauchy tensor σ ij :

σ ij =





σ 11 τ 12 τ 13

τ 21 σ 22 τ 23

τ 31 τ 32 σ 33



 (2.15)

Pressure is always acting along the normal direction relative to the surface of interest and is therefore defined by the diagonal components in the tensor, while the viscous effects act in the surface tangential directions and are identified by the off-diagonal elements. Knowing the stress tensor allows to fully characterize the hydraulic torque, which is generated by the two components. In practical cases it is not possible to characterize every infinitesimal volume with its related stress tensor. This is why it is necessary to find a more suitable way to express the fluid contribution to the overall equilibrium. An alternatives way has been proposed in [7], where the flow effect is treated as a drag force, however, this idea requires to compute a drag coefficient as well, and therefore it does not bring any significant advantage compared to the model by Li and Liou, where the hydraulic torque is seen as a result of a stationary component and a rotational one:

T H = T HS + T HR (2.16)

These two elements can be calculated from their respective coefficients, pre- viously determined by experimental proceedings. The advantage of this formu- lation relies on the fact that through a CFD analysis it is relatively simple to identify the coefficients of interest. The procedure to calculate the static coef- ficient is therefore quite straight forward and will be explained in the following section.

Eq. 2.14 implies also that the weight of the valve is also sensibly affecting the overall force equilibrium. The torque due to gravity can be computed as:

T _w = m[(L × g) · n] (2.17)

where L is the distance between the of gravity and the rotational axis (the

so-called arm length) and m is the disc mass. While L is not varying in its

magnitude, the angle θ formed with a reference plane (from here and afterward

it will always be the vertical one) can change, making the weight torque varying

Figure 2.1: Schematic Representation of the Swing Check Valve

during time. However, the gravity effects do not represent an obstacle to a successful valve motion modeling, since they can be easily defined under every flow condition.

The buoyancy torque is also related to weight and fluid forces; any body immersed in a fluid is subject to a force acting on the opposite direction of the gravity force and directly proportioned to the fluid density. The torque generated by the buoyancy is therefore always opposing the weight torque con- tribution:

T _b = ρ _f V _d [(L × (−g)) · n] (2.18) where V d is the disc volume.

The resulting torque, used also in eq. 1.3 is the submerged disc torque:

T W = T w − T b (2.19)

Friction force is also present in some magnitude, however it can be neglected in particular cases, especially when the contribution is considerably smaller compared to other torque components. This is usually true for big size valves used in industrial applications. In other cases instead, friction factor may not be disregarded, for example in laboratory experiments, when small and light- weight valves are connected to measuring instrumentation. In this thesis work, however, the geometry of interest concerns an industrial check valve designed for nuclear power plant piping system. Its disc weights approximately 60 kg, friction force will be therefore disregarded as it has very little impact compared to weight torque. It is important to note that this simplification will not affect the significance of the results. In [5] the authors themselves mention that fric- tion characterization becomes important only for laboratory tests. A possible formulation for the friction torque is anyway given by [8]:

T f = Z

S

(r × dF f ) · n dS s (2.20)

where dF f is the friction force vector and S s is the contact surface.

Finally, external torques can be applied to the valve in order to maintain it open or close despite the flow conditions, or to reduce closing time. These external torques can be generated by a spring, for example. This torque is usually known in advance and does not present any particular problem.

The right side in eq. 2.14 accounts for the total rotational inertia that is expressed as a sum of different components. In the considered case there are indeed three relevant terms; the first one I is the moment of inertia of the disc itself, generally defined around an axis passing for the a point P :

I P = Z

ρr ² dV (2.21)

Note that ρ can vary as a function of r, which represents the infinitesimal vol- ume dV distance from the axis passing through P. Rotational inertia is usually calculated around the principal axis, i.e. those passing through the center of mass. However it rarely occurs that the rotation is allowed around the center of mass, and in this case it is necessary the use of the parallel axis theorem to find the desired moment of inertia around an axis located at d meters from the gravity center:

I P = I CM + m ass d ² (2.22)

When the disc is swinging, the acceleration of the surrounding fluid should be taken into account in the force equilibrium by the added mass inertial term, I f

in eq. 2.14 . Turesson [7] in his work performed a short test on the importance of including the added inertial component to the equilibrium. The results showed that it is relevant to consider this term in order to achieve a more accurate model. It is important also to note that this added term model is the result of previous works [11] and therefore alternative approaches could lead to different results. For this study’s purpose this added mass model is however adopted, since it is the same used also in [5]. The formula for the inertia of the added mass, which corresponds to the inertia of a sphere having the same diameter as the valve disc and density equal to the that from the fluid is:

I f = 2

5 m sphere R ² _disc + m sphere L ² (2.23) where R is the disc radius and L is the arm length.

Finally, in case of external torques applied to the check valve, the inertia of these components Ie must be included in the equilibrium from eq. 2.14

2.2.1 Calculating Hydraulic Torque Coefficients Accord- ing to Li and Liou’s Model

The Li and Liou model requires to know two coefficients for the stationary and

rotational torque. In a physical experiment presented by the authors, a test

rig capable of extracting the hydraulic torque contribution from the different

monitored moments is presented, and the methodology is explained and tested.

However, this approach is not really convenient, indeed one of the advantages of using a numerical method relies on the relative easiness of finding out the force components generated by the fluid flowing on the valve door. A different approach (from the one proposed in [5]) for calculating the valve coefficients will be applied in this work.

A steady state simulation is firstly needed to identify the stationary torque coefficient C _HS , which will be used afterwards to determine the rotational one C _HR .

The method proposed by Li and Liou for measuring the rotational compo- nent is based on the previous knowledge of all the other torques acting on the disc mentioned in eq. 1.3, such as weight, friction etc.. The equation proposed by the authors to define the hydraulic torque is therefore:

T H = T HS + T HR (2.24)

where:

|T HS | = C HS ρA _v V ²

2 L (2.25)

|T HR | = C HR ρA v

(L ˙ θ) ²

2 L (2.26)

In the previous equations, A _v is the valve disc front area, L is a representative length (the moment arm), V is the reference flow speed, θ is the angular position of the disc and finally ρ is the fluid density.

One of the key points of Li and Liou’s work is dedicated to find a suitable and realizable way for estimating these two coefficients. Several experiments are then performed and described. The basic idea of the researchers consists of performing first a steady state experiment, where the flow speed is constant and the disc position is fixed.The static coefficient is subsequently calculated with equation:

C HS = 2(|T W | + |T E |)

ρA v LV ² (2.27)

Afterwards, the rotational coefficient for transient cases can be computed according to the equation:

C HR = |T W | ± |T E | − (C ˙θ + |T f |) ± |T HS | + (I + I E + I f )¨ θ

0.5ρA v L(L ˙ θ) ² (2.28)

Where (C ˙ θ + |T f |) is the friction contribution and (I + I E + I f ) is the total

inertia of the disc. It is important to note that the procedure illustrated by

Li and Liou to calculate this coefficient is based only on free-falling test with

a steady flow rate. The disc velocity cannot be therefore directly controlled,

but its effect on the coefficient is analyzed by varying it by applying different

external torques.

The CFD simulation presents instead significant advantages, as it is able to compute the hydraulic torque for every case directly from pressure values. There is therefore no need to calculate all the various moment components described in eq. 2.28. Analogously to the physical experiment, a steady state flow simulation is needed to identify the C HS , and in the following transient simulation the stationary torque will be computed using the previously mentioned coefficient.

Another difference between computational and physical experiments relies on the fact that the CFD method is able to not only perform free-falling simulations, but also calculate results for a defined angular velocity of the disc. Both cases will be investigated and compared with the results obtained in [5].

It is also noteworthy that in a further study, an alternative formulation for the rotational coefficient is proposed by Tran [10]:

C HR = C HS cos θ[cos θ + 2V

L ˙ θ ] (2.29)

This expression computes the rotational coefficient value without knowing the hydraulic torque directly for transient cases, which surely is a significant advantage in a physical experiment, as well as in a CFD simulation. However this formulation of the rotational hydraulic coefficient is not the one used for the RELAP5 model, but it would be interesting to compare the results obtained from the two formulations and see if there is a significant correlation. In that case the use of 2.29 should result easier and quicker for calculating the C HR .

Since the CFD code easily computes the hydraulic moment acting on the valve door, the eq. 2.27 and 2.28 become respectively:

C HS = 2|T HS |

ρA _v V ² L (2.30)

C HR = 2(|T H − T HS |)

ρA v L(L ˙ θ) ² (2.31)

In the steady state simulation T HS = T H is directly found by the code, while in the dynamic analysis it is not possible to separate the two parts of the moment [7][8], but T H can still be obtained by the software, while eq. 2.30 can be used to identify the stationary contribution to the overall torque.

2.2.2 Alternative 1D Model for Hydraulic Torque

It is necessary to introduce also an alternative model for the swing check valve dynamics created by Adamkowski [17]. This approach is used in one of the reference works [19] and has been implemented in RELAP5 as well. Some comparison among different models and solvers will be made illustrated later on in the results section.

Similarly to Li and Liou, the Adamkowski model, also referred to as C M −

C _Q , computes the hydraulic torque through some coefficients. The advantage

compared to Li and Liou’s theory is that only steady state experiments are

required. The first of three coefficients express the ratio between measured torque and torque due to static pressure difference across the valve:

C _M = T _H

A s ∆pL (2.32)

The second coefficient, C _Q , characterizes the pressure drop across the valve in a dimensionless way:

C Q = s

ρQ|Q|

2∆pA ² _s (2.33)

The third coefficient C k is simply:

C k = 1

C _Q ² (2.34)

All these parameters are then used in the following equation to compute the hydraulic torque acting on the valve door:

T _H = C _M (θ · sign(Q)) · A _s · ∆P · L (2.35)

∆P = C k

ρQ|Q|

2A ² _s (2.36)

It is important to specify that these coefficients are not unique for every opening stage of the valve, but vary according to the door angle θ. Measuring all the relevant parameters at various angles can therefore describe the valve under a transient flow situation.

2.3 About RELAP5

RELAP5 is a code developed by Idaho National Laboratory and allows one

dimensional modeling of multi-phase flows in piping systems. This code has

been developed initially with the purpose of simulating transients case in Light

Water Reactor systems, but thanks to its generic coding, it can be used to

estimate many hydraulic and thermal transients [18]. The advantage of a 1D

code relies on the fact that complex systems can be rapidly simulated even with

a small computational power, unlike 3D CFD simulations that, for the same

simulation domain, would require much bigger computing resources. This code

presents the possibility of simulating check valves performance, but, in the first

versions of the software, the approach used has been proved not enough accurate,

especially with regard of pressure surges prediction. Some new models have

been therefore implemented, which allows the user to choose among the most

suitable based on the data available. Every model indeed requires specific input

parameters to work, as mentioned in previous sections. This updated version of

RELAP5 that will be used in this thesis work provides 3 choices for swing check

valves; the first one is the original C v model, based on pressure difference across

the valve, which was unreliable. A second approach is the C M − C Q model by

Adamkowski [17]. This method involves the identification of two coefficients

which are then used to compute the hydraulic torque, similarly to the Li and

Liou approach. Finally, the third procedure for swing check valve simulations,

the S r , is based on the Li and Liou C HS and C HR coefficients, using them as

a way to calculate the hydraulic torque. It will be one goal of this work to

investigate the robustness of the newly implemented model in RELAP5, and

compare it to the results given by a 3D numerical experiment performed with

STAR-CCM+.

Chapter 3

Numerical Set-Up

3.1 Geometry

A swing check valve used in a power plant piping system represents the geometry input model for the investigation. In the previous work by Boqvist [8] the same check valve without the lever arm connecting the door to the walls has been investigated, and the analysis results are therefore used as a reference comparison with this work in its early stages. The valve is modeled with a CAD software and then imported into STAR-CCM+, where the geometry of interest is extracted and meshed. The valve is composed by a convergent inlet pipe section with initial diameter of 367.9 mm which ends into a hereby called valve section, that is an approximately spherical chamber where the valve door can swing. This door consists of a disc allowed to rotate around a specified axis (blue component in fig. 3.1). The disc is connected to the valve walls by a lever arm.

After the valve door, the fluid is conveyed again through a contracting section which expands gradually later on and connects to the main pipe again. The allowed valve disc movement is only a rotational one, as mentioned before, and the door position can be completely described by using the angle θ. The vertical plane is the reference: if the disc is perfectly vertical θ = 0 ^◦ . In particular the rotation is constrained by the walls; a complete closed position corresponds to a 4 ^◦ angle θ, and a fully opened one equals to 58 ^◦ . For setting up the simulation properly, it is important to define the flow characteristics. Star-CCM+ has a material database which could be used for standard conditions. The flow characteristic values used in this work are slightly different from the standard settings in STAR-CCM+, and they match those used in [8] displayed in table 3.1

Other than the fluid characteristics, it was mentioned before that moment

of inertia and disc geometry play also an important role when computing the

forces acting on the system. Disc characteristics such as lever arm, mass and

inertia are therefore summarized in table 3.2. These parameters are calculated

with the help of the in-built 3D CAD provided in Star-CCM+

Figure 3.1: Geometry of the Valve

Parameter Value

Type Water

Temperature 176 ^◦ C Op. Pressure 75 bar

Density 883 kg/m ³

Dynamic Viscosity 146 µP a Table 3.1: Fluid Operative Conditions

Parameter Value

Mass 62.2 kg

Density 7859 kg/m ³

Inertia moment around rotational axis 4.841 kg/m ²

Lever arm 255.3 mm

Table 3.2: Valve Door Characteristics

Disc Inertia T ensor =





0.989 0 0

0 0.786 0

0 0 0.483



 [kg/m ² ] (3.1)

3.2 Boundary Conditions

In the original geometry, valve inlet and outlet section are positioned relatively close to the section of interest (i.e. near the valve door). However it is desirable to have a fully developed flow condition for an easier turbulence modeling, this because under these conditions it is possible to define turbulence intensity and length scale in a pipe according to the formulas:

I _turb = 0.16Re ^−1/8 (3.2)

l _scale = D _pipe (3.3)

Extension pipes of approximately 3.5 m on both sides of the valve are suffi- cient to ensure fully developed flow conditions. This extension length is similar to that used in previous investigations. The adequateness of the choice can be easily checked by plotting the velocity profile at several sections in the inlet pipe, this in order to visualize the velocity field and compare it to a typical fully developed profile.

For all the simulations performed in the first part, the boundary conditions applied at valve inlet and outlet are respectively a velocity inlet and a pressure outlet. In particular, velocity inlet sets a specified flow speed and direction for all the cells at the boundary of interest and, for incompressible fluids (like the one considered in this simulations), setting a velocity or a mass flow inlet does not provide any significant change in final results. The outlet is set as a pressure outlet at 0 Pa, this choice allows the possibility to eventually have a reverse flow at outlet, compared to the flow-split outlet type [12].

Other than inlet and outlet conditions, it is important to specify what hap-

pens near wall boundaries. In this case a non-slip condition is applied for all the

walls. Shear stresses are influenced by the surface roughness, and, as mentioned

in section 2.2, they can influence the resulting hydraulic torque. By default

STAR-CCM+ uses the smooth wall conditions to calculate the boundary layer

velocity profile. It is however possible to specify the roughness height for the

wall, which determines the velocity profile, but at the same time a new vari-

able is introduced in the simulation. A sensitivity analysis on the effects of wall

roughness is performed to check if the roughness of typical steel pipes introduces

significant changes compared to the default solution. If the difference is negli-

gible, with the purpose of keeping the model as simple as possible, the default

hydraulic smooth wall is used throughout the simulation. This does not mean

that the boundary layer is disregarded, but it rather implies that the standard

logarithmic-law coefficient is applied for calculating velocity profile across the

layer [12].

Unsteady State Boundary Conditions

In the unsteady simulations there are some differences concerning boundaries conditions. In particular the velocity inlet magnitude will be a function of the simulation time. Other than that, movement constrains are specified for the disc.

More information about the valve degrees of freedom are presented in chapter 3.5.2 . The standard case regards a free falling disc, but other simulations will have a superimposed disc angular velocity expressed as a function of simulation time, similarly to the inlet velocity boundary conditions.

3.3 Convergence Criteria

The CFD software proceeds by iterative steps, which means that the final so- lution is found by refining the results for the governing equations, using the values from the previous step as an input and subsequently correcting them until they converge to the real solution. The procedure allows to get imbalances for every step smaller and smaller. Once the results are not changing anymore in a significant quantity the solution is considered to be acceptable, and has therefore converged. Convergence speed depends on many parameters such as mesh, initial conditions guessing, computational power etc.. STAR-CCM+ dis- plays by default the value of the normalized residuals, which are a useful but not sufficient way to judge simulation convergence. It is important to note that, if the initial conditions for the solution are 100% correct, the residuals will never drop, but the simulation would still be accurate [12]. In the particular case of STAR-CCM+ the software automatically normalise the residuals, which means that for every iteration the imbalance in conservation equations is divided by the maximum value obtained for the imbalance in the first n iterations, where usually n = 5. A different choice of n will therefore give different residual values, and this is why they cannot be used as the only way to judge convergence in a simulation. A common way to decide if the simulation results are acceptable is to monitor a so-called quantity of engineering-interest [12], that could be any physical value relevant for the particular simulation. It is up to the user to define what to monitor during the iteration process. For the steady state simulation in this case it is decided to monitor the C HS in addition to the default residuals.

Another monitored quantity is the total pressure average at the inlet and outlet

section of the pipes connecting to the valve. The coefficient’s stability, iteration

after iteration, is an index of convergence. In the unsteady state simulation it is

harder to check the simulation convergence, as the monitored quantities varies

during time. The technique adopted consist of monitoring the hydraulic torque

acting on the valve door and ensuring that the inner iterations, i.e. the number

of iterative steps performed at every time step, are sufficient to reach a stable

value for every single time step. The stability of the monitored quantity is used

to determine which is the acceptable number of inner iterations to perform at

every time-step. It is important to achieve a good trade-off between accuracy

and time required at every step, otherwise the simulation would take too long

or would be inaccurate.

3.4 Sensitivity Analysis

When performing a computer assisted simulation of any kind, it is important to accurately choose the inputs. The software will indeed provide results for any kind of input, whether it is actually reasonable or not. It is up to the user to critically analyze the results and decide to either accept or refuse them. A typical risky situation is that the final results look nice, everything has converged, the monitored engineering quantities have stabilized, but the values obtained are strongly influenced by the particular settings chosen in the beginning. This can obviously involve misleading results that affects other studies based on them. To minimize this risk, it is a common procedure to perform some sensitivity analysis before running the main simulations. There are several parameters that can be checked for sensitivity and the user has to decide what to check, this because some inputs may be relevant for one case but not for another. Some factors are surely necessary to be tested, while others are expected to not interfere at all in the final results.

It is a common procedure to check the mesh, because different grades of refinement can produce sensible changes in the final results of a given simulation.

In the dedicated section 4.1.1 of this report there is a more detailed explanation on which parameters are tested and why.

3.5 STAR-CCM+ Model Set-Up

The code used in this thesis work is the STAR-CCM+ v. 12.02 developed by CD-Adapco, now part of the Siemens group. The code allows the user to model various physic phenomena such as fluid mechanics, heat transfer, chemical reac- tions. The mathematical models describing the physic continua are formulated on the base of the conservation laws. In this work the phenomenon of inter- est is restricted only to the fluid mechanics and fluid solid interaction which is therefore based on the following equations: conservation of mass, conservation of linear and angular momentum, conservation of energy. The solver converts the expressions to discrete algebraic equations that can be numerically solved, this is done through a discretization method for the continuous domain, either by the finite element or the finite volume method.

3.5.1 Mesh Generation

Static Mesh

In the present work, the mesh generation is performed with the tools provided

by STAR-CCM+, which allows the user to choose among several different kind

of base mesh shapes and refinements. The software allows the choice between

polyhedral, tetrahedral and trimmed cells. In addition to those, it is possible to

introduce prism layers for a better simulation of near wall areas if the boundary conditions allow it . There is also the possibility to use different mesh for dif- ferent regions of the simulation domain, as well as introducing volume extruded meshes in order to extend inlet and outlet boundaries, which is a particularly useful tool for this case.

The procedure adopted by the automated mesher is to first generate a surface mesh, which is basically reproducing the CAD geometry in the best way possible, according to the user defined parameters. Then the volume mesh is created in accordance with the previously surface mesh.

It is really important to have a good mesh quality, that means cells whose skewness angle is acceptable, volume change is moderate, and of course the overall geometry should reflect the original CAD surfaces. There is therefore a built-in possibility to check the mesh validity after its generation, in order to decide whether to refuse or accept the current grid.

Another important aspect to consider when choosing mesh size and type is the time required to solve the equations, that is directly related to the amount of cells composing the grid. A 10 million cells mesh will demand more compu- tational power in order to calculate results within an acceptable elapsed time, compared to a case with coarser grid. This is not only due to the number of cells for which every governing equation has to be solved, but also the convergence rate is affected by that; in a finer mesh the residuals will usually decrease slower [13].

A good mesh is therefore the one that is just enough complex to accurately simulate the physics without using additional unnecessary (and expensive) com- putational power. A sensitivity analysis is performed in order to find out which is the mesh that satisfies this condition. A common way to do that is to create a coarse mesh in the beginning, run the simulation and monitor the desired parameters, then refine the mesh and run everything again. Once the surveyed values are not sensibly changing anymore, the mesh can be accepted and there is no need to use an even finer grid. This procedure has been adopted also for this thesis work, and the results are displayed in 4.1.1.

Except for the cases where it was required to run simulations with different parameters, a base geometry and mesh with defined boundary condition is used for the most cases in the sensitivity analysis. This setup is described in table 3.3:

Mesh Value

Base Size 8 mm Mesh Type Poly Prism Layers 5 Cell Number 2.8 M

Table 3.3: Mesh for Sensitivity and Steady State Analysis

Dynamic Mesh

A different mesh is used in the unsteady simulations, following a method called Overset Meshing, which is really useful for simulating the closing valve disc.

However, this methodology requires a slightly different mesh setup, and an additional sensitivity analysis is performed as well.

The Overset mesh consists always of at least two regions partially overlap- ping each other. The proper definition of the boundaries and interfaces allows the solver to divide the cells into active, inactive and acceptor cells [12]. As the name suggest, the inactive cells are deactivated and the flow is not solved for them while the acceptor cells are used to couple the solution between overset and background region. This meshing technique allows not only to simulate a moving body with up to six degrees of freedom, but it also permits the creation of a zero-gap interface, which is of particular interest when simulating the full closure of the valve, avoiding simulation crash when the geometry does not allow the fluid to flow anymore. This was an issue in previous works, it was almost impossible to compute the fluid properties with the valve closed. With the zero- gap wall, the code treats a specified amount of cell layers as a wall boundary when the gap becomes smaller enough. In this way, the flow sees a virtual wall when the disc gets very close to the seat. Other than this, the overset mesh allows the disc movement thanks to the interpolation of value between overset and background region. As the name suggests, the overset region is allowed to move into the background region, and the solver recognizes automatically when to use either the overset or the background cells for the solution. The cells in Overset region that are overlapping the background domain will be treated as active or acceptor, while the external cells will be deactivated. In this way, if the Overset region moves, the cell status will be updated and eventually some cells will become active while other will be excluded from the solution domain.

This is exactly the case for the closing valve disc, which is part of the overset region, rotating according to the momentum exerted by all the forces described in the theory section. While the disc is moving, some cells in the background becomes inactive because the door is occupying that space. Later on the same cells will become active again since the disc has changed position.

The main disadvantage of using this technique relies on the fact that an additional interpolation has to be performed at every iteration, as well as the fact that the moving disc requires the update of the Overset-Background interface for every time step. This all results in increased time required for the computation, especially because the Overset algorithm is not fully optimized for parallel solver processes, which means that even if the solution is computed on several cores simultaneously to speed up the procedure, the benefits will not be so significant as they could.

The importance of choosing a mesh which is not too fine is therefore even

more relevant in this case, and a proper sensitivity study is performed to obtain

the best solution in terms of accuracy and solving speed.

Background Overset

Base Size 8 mm 7 mm

Mesh Type Poly Poly

Prism Layers 5 5

Cell Number 2.65 M 0.56 M Table 3.4: Final Mesh for Unsteady State Analysis

3.5.2 Fluid - Body Interaction

As expressed in detail in section 2.2, the valve disc movement is governed by Newton’s second law for momentum. The STAR-CCM+ solver is able to com- pute all the forces acting on the body, in particular those due to the fluid inter- action and the weight. This built-in solver function is called 6-DOF (Degrees of Freedom), and can fully describe the movement of a rigid body interacting with the flow. Additional external forces can be also defined in order to simulate for example springs connected to the body, or any sort of damping system. By knowing the inertia tensor, the initial center of gravity position and the body mass, the 6-DOF solver can find the position, velocity and acceleration of the body at any time step.

In this particular simulation case, the body has only one degree of freedom, which is the rotation around the pivoting axis. Therefore the only quantity com- puted by the 6-DOF solver is the total momentum acting around the rotational axis. In addition to that, it was mentioned that the disc cannot open more than 58 ^◦ , and cannot be less than 4 ^◦ either. These limitations are implemented in the solver, which will add a fictitious momentum to avoid further movement when the disc reaches the extremes of the allowed angles interval.

The use of this solver not only allows to describe the body motion itself, but it is also provides necessary quantities when computing the C HR , as it is required to know the rotational velocity of the disc, according to eq. 2.31.

3.6 RELAP5 Model Set-Up

3.6.1 Check Valve Modeling Input

As briefly introduced previously, one of the main goals of this thesis work is to compare the robustness of a one-dimensional model compared to a 3D numeri- cal experiment. The mono-dimensional calculations are performed by RELAP5.

This software accepts a specifically formatted input file which represents the modeled system and provides also the initial conditions. In the recently im- proved version of RELAP5, the swing check valve can be modeled according to [5] and [17], provided that the required input coefficients are available. These can be supplied by the user in form of control variables i.e. user-defined quan- tities, or as tabulated values. [15]

Control variables in RELAP5 represent a way to define single variable func-

tions and algebraic expressions directly into the input code. This approach, in this specific problem formulation, allows the user to define two functions which describe a coefficient each, making in this way possible to implement the regression results explained in sections 4.1.2 and 4.1.4.

The tabulated-values way presents instead a slightly different procedure to use the CFD analysis results. While using a recently implemented feature called freetab, the software automatically performs a linear or polynomial interpolation between the user-provided data. However, the table has to be formatted in a specific way, which is not always the best or easiest solution. Considering for example the rotational C _HR coefficient, the input table needs indeed a series of three vectors, one for each independent variable θ, ω, |V |, plus a fourth vector containing the coefficient values, whose length is given by:

length(C HR ) =

3

Y

i=1

length(V i ) (3.4)

where V i indicates one and each of the three vectors defining the independent variables. This means that a table can easily contains several thousand elements, and has to be therefore provided as an external file, since the in-built table format can accept up to 99 entries. If the freetab approach is chosen, it is necessary to use an interpolating function to create the required table. This function can be obtained by CFD results, as discussed in the proper section in 4.1.4.

In the model implementation for this thesis work, however, it was chosen to proceed by using the regression functions defined for different opening angles, explained as well in 4.1.4. This because the large amount of data available from CFD simulations made it possible to create a good model. In addition to that, the 3D interpolating function in MATLAB did not provide consistent values for certain combinations of independent variables, meaning that the freetab solution presented unrealistic results for some valve angles and mass flows combinations.

This is probably due to the difficulties of fitting the data collected while the disc was almost completely opened or shut. At those points indeed the computed C HR assumed inconsistent values due to, for example, the extremely low angular speed, or the extremely high flow velocities around the disc when it was about to reach the valve seat.

The scattered interpolation forced the function to pass through every single collected data, which means that every point had the same importance in the final result. With the generalized regression model instead, the points that are laying far away from the general trend are still influencing the overall model, but their weight is reduced, giving a more reliable trend. This is the main reason why it was chosen to proceed with this alternative approach; the implemented code used in RELAP5 can be found in appendix B. Further details about the regression model are explained in 4.1.4.

Besides the coefficient values, RELAP5 requires additional information to

be able to properly characterize a 3D object in a 1D way. In particular, it is

necessary to specify some relevant angles to perform the calculation (β and θ in

figure 3.2). In addition to that, it is also necessary to specify other parameter

Figure 3.2: Angle Reference System

such as valve door volume, lever arm length (R CG ), disc weight etc.

To achieve the best resemblance between CFD and RELAP5 simulation is it therefore extremely important to correctly set all these parameters, allowing the code to compute an accurate solution.

3.6.2 System Description

The first system simulated in RELAP5 consists of pipeline composed by a simple

set of pipes connecting a pump and the check valve. The pressure in the system

is defined at the inlet and outlet, the mass flow is therefore determined according

to the pump curve details defined in the input. The basic RELAP5 code used

for the simulation was obtained by a previous thesis work by Eriksson [19], but

in that case the check valve was modeled according to [17]. After some changes

in the input file (see appendix B), it was possible to use the new model [5]. In

the simulation it is possible to monitor what happens when the pump is shut

down, in particular it is interesting to check the mass flow through the valve,

and the closing time required. The system here presented is however witnessing

a relatively slow transient, indeed, the valve closure happens in approximately

0.15 s, with the flow acceleration that is −11 m/s ² on average. However, in some

situations like a pipe ropture, the flow velocity can change much faster, with

accelerations reaching up to −400 m/s ² . It is therefore interesting to check if

RELAP5 can properly handle even these extreme situations. To do that, a new

simple pipeline where the user is able to control the mass flow is defined. The

system analyzed consists of a series of long pipes, a time dependent junction,

which ensure the possibility to directly control the mass flow rate acting as

piston pump, and the swing check valve.

Figure 3.3: Scheme for the First Pipeline

Figure 3.4: Scheme for the Second Pipeline

Chapter 4

Results

Both sensitivity and sharp analysis are performed during this thesis work, and this led to a consistent amount of data to present and re-elaborate. The use of different softwares for 1D and 3D analysis requires a structured presentation of the results. First, the CFD analysis of the swing check valve is presented, where the most significant results surely regards C HS and C HR characterization.

Then, the RELAP5 simulations are presented and finally again the 3D CFD comparison results with the previously mentioned code are displayed.

4.1 Valve Investigation with 3D CFD

In this first part of the work, the focus is posed on the valve characterization with reference to [5]. In particular the performed analysis aims to identify the two coefficients according to the Li and Liou’s model introduced in chapter 2.2, as well as the three other coefficients for the Adamkowski’s model.

4.1.1 Sensitivity Analysis For Steady State Simulations

Sensitivity analysis is essentially the procedure of monitoring some interesting parameters, which can help the user decide for a model rather than another. In order to present a clear and structured analysis, in every subsection hereafter treated, results show the C HS variation and the pressure drop coefficient ψ.

Besides meshing techniques, several other factors have to be defined when

setting up the simulation, and the influence of these on the final result cannot

be known in advance. As mentioned before, the ultimate goal is to achieve a

simulation setup which provides a good physic resemblance of a real experiment,

but at the same time it is desirable to keep it as simple as possible. This

because every additional input affects the output, and an excessively complex

model has more chances to encounter errors. Several analysis are therefore

performed also on other settings, in particular, studies are conducted on wall

roughness, turbulence models, unsteady vs. steady state simulation, Reynolds