T Thénint

E Balmes
email: balmes@sdtools.com

M Corus

STABILIZATION EFFECT OF SHOCK NON-LINEARITY ON THE DYNAMICS OF A STEAM GENERATOR TUBE

Keywords: dynamics, shock, stability

In Pressurized Water Reactors of nuclear power plants, steam generators act as heat exchangers between primary and secondary coolant fluids. They consist of a bundle of U-tubes in which flows the primary coolant fluid. Several support plates guide these tubes vertically. Secondary coolant fluid flows along the U-tubes and passes through the space between tubes and plates. This space, initially of a foliate shape, is filled with sludge deposits. Consequently, fluid flow is accelerated and the tubes are more excited. Moreover, the mechanical bonding between tube and plate is changed. The combination of these two phenomena can lead to dynamic instabilities and tube cracks. This paper focuses on the effects of shocks on tube instability. As we concentrate on a mechanical point of view, the study of the nonlinear dynamics is made in the following specific conditions: the tube is in air with no flowing fluid; an instability is generated by injecting a force proportional to the velocity of a point of the leg of the tube. The tube is modeled as an Euler-Bernoulli beam. The contact between tube and plate is computed assuming circular obstacles and a contact force that is linear with a gap. A reduced model is generated to represent the bandwidth of interest and the effect of contact forces. Time evolution is then computed using a nonlinear Newmark scheme. Numerical simulations show the effect of shock nonlinearities on the dynamics. Cases of unstable unconstrained tubes leading to bounded stable dynamics when shocks occur are analyzed. The instability in the unconstrained condition does not imply that the bilateral contact condition is unstable. Areas of the contact stiffness/feedback gain plane are shown to lead to stable modes of the bilateral model. The work then presents results on an experimental bench that includes the curved part of a tube and one Support Plate. The instability is generated by feeding back to a shaker the amplified signal of a velocity measurement performed with an integrated accelerometer. For great enough feedback gains, the vibration amplitudes grow until shocks occur. The permanent regime eventually found is analyzed. A good similarity between test and analysis is found.

INTRODUCTION

In Pressurized Water Reactor (PWR) of nuclear power plants, steam generators exchange heat between primary coolant fluid and secondary coolant fluid (Figure 1). They consist of a bundle of inversed U-tubes in which flows the primary coolant fluid. These tubes are clamped at the bottom of their legs. Several support plates maintain them vertically. Secondary coolant fluid flows along the tubes and through the holes made on the support plates. The space between the tubes and the support plates, initially of a foliate shape, is filled with sludge deposits. Consequently, fluid flow is accelerated and the tubes are more excited. Moreover, the mechanical bonding between tube and plate is changed, from free gap to clamping. The combination of these two phenomena can lead to dynamic instabilities and tube cracks.

Different methods are available to predict tube instability in regard with external fluid flow. The method developed by EDF uses linear boundary conditions to model the interaction between the tube and the support plates, typically clamping or supporting. Investigations have been made to get better understanding of the role played by shock non-linearity in respect with tube instability.

Piteau et al. [START_REF] Piteau | Vibro-impact experiments and computations of a gap-supported tube subjected to single-phase fluid-elastic coupling forces[END_REF] have made experiments and computations of a vibro-impacting straight tube subjected to fluid-elastic forces and planar displacements. Shocks modify the instantaneous response frequency and thus modify the fluid-elastic forces, explaining how shocks can stabilize this unstable dynamics. Here, we investigate the stabilization effect of shock nonlinearity by itself. The study of the nonlinear dynamics is made in the following conditions: the U-tube is in air with no flowing fluid; instability is generated by injecting a force proportional to the velocity of a point of the leg of the tube.

Section 2 presents the experimental setup. Contact model and reduced basis used for numerical simulations are presented in section 3. The stability of a tube with linear contact stiffness is studied in section 4. Finally, numerical and experimental non-linear simulations are presented in section 5.

EXPERIMENTAL BENCH

The experimental bench consists of a portion of a U-tube (Figure 2), constrained at its basis to be as near as possible to a clamped condition. The two legs of the tube are inserted in fo-liated obstacles, whose geometry reproduces real support plates. These obstacles have been machined with different radius, the gap between the plate and the tube approximating the quantity of sludge deposits. The maximum gap available is around one percent of the external diameter of the tube. The thickness of a support plate is about 1.5 times the external diameter of the tube. A shaker excites the leg in the out-of-plane x direction (z is the vertical and y is an in-plane direction). The shaker is fixed on a lift table, so as to adjust the height of the destabilizing force. Twenty-two B&K 4375 sensors measure accelerations at different points of the structure. A B&K 2635 charge amplifier integrates the acceleration and feeds an electrical voltage proportional to the velocity of the excitation point. Measurements are sampled using LMS Scadas III Front-end, driven by LMS Test.Lab software.

NUMERICAL MODEL

A reduced model has been built to perform non-linear numerical simulations. An updated numerical model of the experimental bench is created in section 3.1. The contact model between tube and plate is presented in section 3.2. Reduced basis is built in section 3.3. Numerical simulations are compared with test results in section 3.4.

Beam model

The tube is modeled as a curved beam. The length of the elements is compatible with a maximum frequency of interest linked to the spectrum of the contact loads. A preliminary study of the kinematic in the contact area has led us to consider an Euler-Bernoulli formulation rather than a Timoshenko formulation.

The mechanical connection between the tube and the frame is quite stiff but is far from clamped conditions. We represent this bind with six discrete springs, three affecting the trans-lation degrees of freedom and three affecting the rotation degrees of freedom. A modal analysis gave us twelve modes under 200Hz. The tube material parameters and the stiffness of the discrete springs were updated taking into account the frequency of the modes 4, 7, 9 and 11, out-of-plane bending or torsion modes.

The frequencies smaller than 100Hz are quite well correlated (Table 1). The MAC (Modal Assurance Criterion), normalized scalar product between two vectors, is very close to 1, indicating a linear relation between experimental mode shapes and updated numerical mode shapes. Mode shapes are quite close to the ones of a base-clamped tube (Figure 3). A double mode around 85Hz has been badly identified experimentally and the updated numerical frequencies and deformed shapes are not consistent.

Mode number

Measured As the mode frequencies are well separated and providing that the damping is small, we make the modal damping assumption. The few well identified damping values will be used, typically 1.5% for the two first modes and 0.4% for the third mode. The other modes will be affected with 0.2%, common value for the metallic structures cast in one piece. The high frequency modes above 7500Hz (this limit is derived from measured contact forces) will be arbitrary affected of a 2% damping ratio.

Contact model

The contact between tube and plates is solved using a functional representation of contact forces: if the displacement is smaller than a defined gap, then the contact force is null; else, the contact force applied to the node is proportional to the penetration and normal to the sur-face of the obstacle. The global stiffness of the contact between tube and plate has been experimentally measured at 4.10 6 N/m, which is of the same order than the ovalization stiffness of a pipe presented in [START_REF] Antunes | Méthodes d'analyse en dynamique non-linéaire des structures[END_REF] or than the stiffness used in [START_REF] Piteau | Vibro-impact experiments and computations of a gap-supported tube subjected to single-phase fluid-elastic coupling forces[END_REF].

The geometry of the obstacle has been simplified, assuming that circles would be a good approximation of the foliated holes. The stiffness has to be distributed over the height of the support plate only using a few planar obstacles (Figure 4). This is done degrading the initial contact where both the plate and the tube are modeled as solids. The contact between a beam and the plate would be modeled as the interaction between a neutral fiber and a contact stiffness density. Here, contact conditions are assigned to nodes interacting with discrete springs whose value is the global stiffness divided by the number of planar obstacles. When the gap between tube and plate reduces to zero, the translations along x and y are blocked. Due to the support plate thickness, the corresponding rotations are also blocked. To block the rotations, at least two levels of obstacles must be considered. When the gap is zero, the tube is connected to the discrete springs. The evolution of the eigen modes with respect to the contact stiffness is compared for different number of discrete obstacles. There are twentyone nodes along the plate thickness so at most twenty-one obstacles can be considered. Five levels of obstacles imply very limited errors on frequencies compared to those found with twenty-one obstacles. No shock damping was used as no physical value has been identified. Friction was not taken into account. We consider that the dissipation induced by the modal damping stands for all the physical dissipation.

Reduced basis

The tube non-linear dynamics are computed solving the discretized equation:

shaker shock ( ) ( ) NL MX CX KX F t F X + + = + && & (1)
where M, C and K are the mass, damping and stiffness matrices of the unconstrained tube. The external excitation field has been separated between the excitation term that stems from the shaker and the contact/impact forces located at five nodes per support plate. To reduce the size of this system, we project (1) on a reduced basis T. The degrees of freedom q are solution of (2), with X Tq = :

shaker shock shaker shock ( ) ( ) ( ) T T NL NL r r r M q C q K q T F T F cTq f t f q + + = + = + && & (2)
and [c] is an observation matrix used to express the penetration as a linear combination of the degrees of freedom. We now detail the reduction basis building procedure. The basis must be rich enough to solve precisely the non-linear dynamics but should be as reduced as possible to allow quick computations. The basis must describe tube displacements when the gap is large, i.e. the tube is unconstrained, and when the gap is close to zero, i.e. when the tube is constrained. Moreover, shock forces applied at a few nodes must be well projected. Finally, the induced spectral bandwidth must include the high frequencies excited by contact/impact forces (around 7500Hz).

The quality of the reduced basis is first based on the computation of the eigen modes of the tube with bilateral contact conditions. We compare frequencies and deformed shapes found solving the reduced problem or solving the un-projected problem. The error is not significant for a wide range of bilateral contact stiffness if the reduced basis is built concatenating deformed shapes of:

-the unconstrained tube, -the tube with x and y translation degrees of freedom blocked over the thickness of the support plate, -the tube with x and y translation degrees of freedom blocked at middle node of the support plate, -static corrections which correspond to unitary force imposed on each degree of freedom of the obstacles (KX=b i ). Family of vectors T is no more orthogonal with respect to M and K. That is, some vectors are nearly collinear. For example, the static corrections are nearly collinear one to each other and nearly collinear to the first in-plane and out-of-plane bending modes. Collinear vectors are eliminated and the remaining family is orthonormalized with respect to M and K, as presented for example in [START_REF] Sternchüss | Multi-level parametric reduced models of rotating bladed disk assemblies[END_REF]. The first deformed shapes, corresponding to the unconstrained tube, are not changed. In fact, only the information (frequency and deformed shapes) about high frequencies has been modified. As C r is assumed diagonal by hypothesis (modal damping), all the matrices appearing in (2) are diagonal. Thus, the numerical computations are accelerated.

The choice of the reduced basis, presented in Table 2, has been validated with non-linear numerical simulations. "Blocked modes" and static corrections are really useful: the same quality of results is obtained but using fewer modes.

Comparison of numerical simulations with test results

Numerical simulations have been computed for different gaps and different imposed excitations, without velocity feedback. A non-linear Newmark average acceleration scheme [START_REF] Gérardin | Théorie des vibrations. Application à la dynamique des structures[END_REF] is used, with a time-step of 2.10 -5 s, small enough to avoid energy errors during non-linear time integration. No numerical dissipation is introduced. Measurements are performed at a sampling frequency of 2048Hz. Power Spectral Densities (PSD) (see [START_REF] Kenneth | Vibration testing: theory and practice[END_REF] for example), are computed. The PSD estimator is normalized to compare signals with different time-steps: A random voltage filtered between 10 and 40Hz is fed to the shaker; the injected force is modulated by the interaction between the tube and the shaker. This force is measured and used as excitation in the numerical simulations. Accelerations are computed or measured at the tube apex, in the out-of-plane direction (Figure 6). Test-analysis correlation is correct over the frequency band [10Hz, 200Hz]. For most peaks and most levels of excitation, amplitudes and frequencies match. The spread of the peaks around 90Hz and 170 Hz is well represented. The anti-resonance around 120Hz is both measured and computed. But the emerging peaks around 28-30Hz and 70Hz are not well reproduced numerically, with the response overestimated and frequencies that not really match.

The satisfactory correlation between test results and numerical simulations validate the construction of the numerical model. A beam model, a contact model, a reduced basis and an integration scheme have been chosen to accurately and quickly compute the non-linear dynamics of a U-tube impacting support plates.

EFFECT OF BILATERAL CONTACT ON COMPLEX POLES

An apparent contact stiffness can be derived from the retained shock model by dividing the contact force with the node displacement. This apparent contact stiffness is null when the contact force is null and tends to the contact stiffness when the force goes to infinity. At every time instant, depending on contact conditions, an underlying linear system exists. In this part, we consider a very simple underlying linear system: each contact node is affected bilateral contact conditions with equal contact stiffness, varying from 0 to 10 8 N/m, as if the gap was null.

Assuming that a feedback loop makes the shaker inject a collocated force proportional to the out-of-plane velocity of a tube leg point and replacing contact/impact forces with bilateral contact conditions, the dynamics ( 2) is now solution of:

[ ] [ ] T T r r r bil loop bil M q C q K q g T c c Tq k K q g C q k K q + + = - = - && & & & (4)
with [c] the observation matrix of the degree of freedom whose velocity is measured and g the gain of the feedback loop; [K bil ] is the unitary bilateral contact stiffness matrix and k is the apparent contact stiffness. It is worth noting that C loop and K bil are not diagonal matrices. We study the stability of this linear system with respect to the position [c], gain g and stiffness k.

To do this, we compute the complex eigenmodes of the following system:

( ) ( ) { } { } 2 0 i r i r loop r bil i M C g C K k K λ λ ψ   + - + + =   (5)
The complex poles are expressed in terms of frequency ω i and damping ratio ζ i (6). They are arbitrary sorted by increasing frequency, which is a convenient but delusive choice.

2 1 i i i i i j λ ζ ω ω ζ = - + - (6)
We are interested in destabilizing low frequency modes. This is done choosing the point where to inject the feedback force, about 20cm below the support plate. The modes which can be easily observed, i.e. with high [cT i ], are the ones whose damping ratio is changed. Every mode with large out-of-plane components becomes unstable except for the second torsion mode and the third out-of-plane bending mode (Figure 7): a vibration node is located near the shaker. We remark that frequencies are not modified when g increases. Damping ratio is a linear function of g. Even if the damping ratio is quite high (-10% for the first mode and g=20), mode shapes are very correlated to the mode shapes of the unconstrained model. Pole evolution in the complex plane with respect to the apparent contact stiffness is much more complex, with many crossings. In the meantime, deformed shapes are continuously transformed. For large enough apparent contact stiffness, the deformed shapes are no longer correlated to those of the unconstrained model. The MAC matrix between these two sets of deformed shapes is not diagonal (Figure 8). We focus on the evolution of frequencies and damping ratios (Figure 8), for a given feedback parameter (here, g=5). Frequencies increase with respect to the apparent contact stiffness. Damping ratio evolution is more complex. For some modes, damping ratio is not affected by the apparent contact stiffness. These modes were stable and stay stable. The first in-plane flexion mode has an initial damping ratio of 1% which decreases to 0.2%. Initially unstable modes become stable, their damping ratio increasing with respect to the apparent contact stiffness. The mechanical system is stabilized for apparent contact stiffness greater than 10 6 N/m: every damping ratio is positive.

a) b)

Figure 8 Evolution of the frequency (a) and the damping ratio (b) with respect to the apparent contact stiffness. Green points denote stable systems, with no negative damping ratio.

For a given feedback parameter, adding bilateral contact conditions leads to pole modifications and unstable modes are transformed to stable ones. The stability of the mechanical system has been studied for a wide range of feedback parameter and apparent contact stiffness (Figure 9). If g is small enough, the system is stable. If g is large enough, the system is unstable. If not, the system is unstable unless k is greater than a given value, which increases as g increases.

Figure 9 Stable configurations with respect to the feedback parameter and the apparent contact stiffness.

NON-LINEAR SIMULATIONS

Adding bilateral contact conditions to the destabilized tube can lead to stable mechanical system. In this part, we investigate the effect of shocks on the stability of the tube, subjected to a positive feedback loop on the velocity. Numerical and experimental results are presented. Computed out-of-plane displacements and measured accelerations are observed (Figure 10). Support plates are between P3 and P4 and between P8 and P9, i.e. close to the top of the legs. The shaker is right under P3.

Numerical results

Initial conditions are null. Tube is excited at its apex by an impact or by a low-pass filtered noise. Non-linear simulations are made over 30s, with a time-step of 2.10 -5 s.

For low values of g, displacements are bounded (Figure 11). For g=0, free vibrations of a damped mechanical system are observed. g=6 should imply an unstable dynamics, as the apparent contact stiffness stabilizing it is greater than the contact stiffness used in the contact model. Feedback loop brings energy to the tube, which is dissipated during contact phases and in high-frequency modes. The vibration regime is stationary. When displacements grow, penetrations and contact forces grow. This can define a criterion describing the stability of a computed non-linear dynamics: if contact forces exceed 1000N, which corresponds to a penetration equal to the tube thickness, the simulation is described as unstable. Non-linear contact forces computed in section (3.4) were less than 50N. For greater values of g, displacement amplitude increases (Figure 12). Even if crest-to-crest out-of-plane apex displacement reaches 8mm, which is quite high, simulated dynamics with g=12, seems stable. When displacements are big enough, contact forces grow and energy is dissipated. Unstable dynamics is exhibited for g=12.1. Contact forces do not stabilize the system, damping stays negative and amplitude of displacements grows until maximum contact force criterion is met. The same behavior has been observed when the impact excitation was replaced by a lowpass filtered noise. Feedback parameters lower than g=12.6 lead to stable dynamics. Contact forces can numerically stabilize unstable dynamics, and for a wider range of feedback parameter g than a bilateral contact model.

The apparent contact stiffness distribution drives the maximal stable value of g. For nominal value of gap between tube and support plates, bilateral contact model with equal stiffness affected at each discrete obstacle is a bad approximation of the real distribution of apparent contact stiffness (Figure 13). Apparent contact stiffness is not uniformly distributed over the different obstacles. Superior obstacles apparent stiffness is globally higher than this of inferior obstacles. The apparent stiffness of the left leg obstacles is totally linearly uncorrelated to apparent stiffness of the right leg obstacles. The apparent stiffness of neighbour contact points is highly correlated and decreases with the distance.

Experimental results

The destabilization is achieved by feeding to the shaker a velocity signal obtained by integrating measured acceleration. The feedback parameter g is fixed using a voltage amplifier. Due to the coupling between tube and shaker, the injected force is not simply proportional to the velocity: peak amplifications and differences of phase are observed. To impose null initial conditions, shaker rod and tube are manually stopped and then let free. Measurement noise is thus amplified by the feedback loop and excites the tube. If g is large enough, vibration amplitude grows. Shocks occur and vibration amplitude remains bounded.

Different regimes can be distinguished in Figure 14. The initial amplitude growth clearly appears on all measurement channels. A stabilized shock regime lasted for about one minute. Then, three different regimes where observed. They differ by their spectral repartition of energy, peak frequencies not being modified.

Operating Deflection Shapes (ODS) can be exhibited. They have vibration nodes located at nodes close to the support plates. They differ from mode shapes of the unconstrained model. They look like mode shapes of the bilateral contact model: ODS are combinations of different shapes, with phase difference, for example a torsion mode combined with an in-plane bending mode. In-plane movements are of the same amplitude as the out-of-plane movements: the specific shape of the foliate hole couples these directions. Despite drawbacks on velocity measurement and on injection of a proportional force, the tube has been destabilized, and it has been shown that contact forces kept vibration amplitude bounded.

CONCLUSIONS

This paper presents the stabilization effect of shocks on an unstable mechanical system. An experimental bench consists of a portion of steam generator U-tube. A feedback loop injects a force proportional to a measured velocity. The tube is destabilized by this "negative damping". Foliate obstacles stop the growth of vibration amplitudes, which remain bounded. A numerical model has been built. The choice of a beam model, a simplified contact model and a reduced basis has been presented. Numerical simulations are quick and correlated to measurements, in the absence of feedback loop. Feedback loop is idealized for numerical simulations: injected force is proportional to the velocity. Unstable free systems are kept bounded by shocks. A possible explanation of the stabilization effect is given studying a linear bilateral contact system: contact stiffness modifies the complex poles. Some work has still to be done to make the three approaches converge and predict the same feedback parameter leading to global instability and dynamical divergence. The contact model should be more complex, the specific foliate shape coupling in in-plane and out-of-plane motions. The feedback loop from measured velocity to injected force should be modeled and affected to numerical simulations. Finally, bilateral contact is a too rough linear model: apparent contact stiffness is not uniformly distributed over the obstacles.

  Figure 1 a) Scheme of a PWR nuclear plant. Three fluid circuits are described: primary circuit (red), secondary circuit (blue) and a cooling circuit (green). b) Scheme of a steam generator. The bundle of Utubes and support plates are shown.

  Figure 2 a) Global view of the experimental bench. We distinguish the two foliated obstacles right under the curved part and the shaker on the right leg. b) Detailed view of a foliated obstacle.

3

 3 Figure 3 First four modal shapes of the updated numerical model: out-of-plane bending (a), in-plane bending (b), torsion (c), out-of-plane bending (d).

4

 4 Figure 4Different contact models: a) solid-solid, b) neutral fiber-stiffness density, c) beam nodesdiscrete springs.

Figure 5

 5 Figure 5Evolution of the first eigen frequency with respect to the global contact stiffness. The modes go from unconstrained boundary conditions to supported and then clamped ones.

M

  Figure 6 PSD of the out-of-plane apex acceleration measured (green) or simulated (blue), for increasing RMS value of the excitation (a to d). Vertical lines indicate unconstrained modes of the tube.

  Figure 7 a) Evolution of the complex poles with respect to the feedback parameter. If the real part is positive, the mode is unstable. b) MAC-M between the unconstrained modes and the complex modes obtained for g=20.

  Figure 10 a) Sensor location. b) Corresponding out-of-plane translation degrees of freedom. For example, P4X corresponds to 471.01. The shaker is indicated with a red point.

  Figure 11 Computed displacements (blue: 236.01; green: 261.01; red: 471.01; cyan: 426.01) for different values of g. a) Damped vibrations. b) Stationary vibrations: displacements are bounded.

  Figure 12 Computed displacements (blue: 236.01; green: 261.01; red: 471.01; cyan: 426.01) for different values of g. Tube is initially excited by an impact. a) displacements grow but remain bounded and eventually decrease; b) unstable regime: contact forces do not stabilize growing displacements

  Figure 13 a) Cumulative distribution function of apparent contact stiffness. Contact nodes from bottom to top of the support plate: N153 and N438 (blue), N158 and N443 (green), N163 and N448 (red), N168 and N453 (cyan), N173 and N458 (magenta). b) Statistical correlation of apparent contact stiffness.

  Figure 14 a) Measured acceleration signals. Red lines separate the different observed regimes. b) PSD of the stabilized regime (16s to 78s). c) ODS at 37.8 Hz for the stabilized regime: torsion mode. d) ODS at 230.7 Hz for the stabilized regime: out-of-plane bending.

Table 1

 1 Comparison between measured and updated frequencies. Damping ratio and MAC.

		frequency	Measured damping	Updated frequency	MAC
		(Hz)	ratio (%)	(Hz)
	1	4.9	1.5	4.9	0.97
	2	8.5	1.5	8.4	0.97
	3	16.5	0.4	16.7	0.99
	4	31.8	0.2	31.7	0.98
	5	32.2	0.2	32.1	0.97
	6	47.8	0.2	47.5	0.99
	7	49.9	0.2	49.6	0.98
	8	88.5	0.2	87.7	0.80
	9	88.9	0.2	88.4	0.66

Table 2

 2 Comparison of different projection basis. Orthonormalisation has eliminated half of the vectors and affected high frequencies to static corrections and blocked modes.

ACKNOWLEGEMENTS

The authors want to thank I. Negreanu (EDF R&D) for his valuable contribution on experimental setup and tests.