Ramiro Mena

Stéphane Guinard

José V Aguado

Antonio Huerta

Distortions in Large Aluminum Forgings: Current Situation and Simulation Challenges

Keywords: Distortions, Residual Stresses, Reshaping, Simulation

Distortions after machining of large aluminum forgings are a recurrent problem for the aeronautical industry. These deviations from design geometry are caused by the presence of residual stresses, which are developed along the manufacturing chain. To solve this problem, a series of post-machining operations called reshaping are required. Despite reshaping manages to restore the correct geometry, it is highly manual and time-consuming, therefore, there is a need at an industrial level to use numerical simulation to study reshaping. The present document describes the problem of distortion, the operations required to mitigate these geometrical defects and the challenges associated to simulate reshaping.

Introduction

Large and thick-walled aluminum forgings are widely used in the aeronautical industry. The key material properties as a great strength-to-weight ratio and good formability in combination with the manufacturing process allow producing complex shapes in an economical way [START_REF] Kalpakjian | Manufacturing Processes for Engineering Materials[END_REF]. Additionally, when compared to other metal working processes (e.g extensive machining, welding or casting), improved material properties as grain size and orientation are obtained [START_REF] Prasad | Processing of Aerospace Metals and Alloys: Part 2 -Secondary Processing[END_REF]. However, aluminum forgings show important distortions after machining due to the presence of residual stresses which are a consequence of the non-homogeneous plastic strains developed along the manufacturing chain, especially after the heat treatment of quenching [START_REF] Ma | Modeling of Residual Stress and Machining Distortion in Aerospace Components[END_REF][START_REF] Robinson | Residual stress in 7449 aluminium alloy forgings[END_REF].

Such distortions are variable in nature, but the industrial process incorporates a robust post-machining stage called reshaping. Particularly, skilled boilermakers operate a sequence of mechanical loadings on the non-conforming part and re-establish its nominal geometry while obeying guidelines imposed by stress offices -which guarantees the produced part has still mechanical characteristics compatible with its operational role -. An example of the structural parts subjected to reshaping are the cruciform beams, located at the wing-box in an aircraft and highlighted in red in Figure 1. On the other hand, reshaping is a highly manual and collaborative process as multiple operators are involved to repair one single part. At the same time, the process fully depends on the experience of the boilermakers and repairing each structural part in a tailor-made manner is time-consuming. Therefore, there is a strong demand from the industry to derive operational assistance from capabilities offered by numerical simulation tools [START_REF] Sim | Challenges of residual stress and part distortion in the civil airframe industry[END_REF].

Focus on distortion mitigation via reshaping, numerical simulations may offer huge margins for improvements such as operational assistance to boilermakers and related gains in time-cycle (optimized reshaping sequence, optimized parameters for each elementary reshaping operation), improved knowledge of actual part and related gains in safety (monitoring of actual plastic strains). However, some research efforts are required in order to migrate from the actual heuristic repairing method to a simulation-assisted operation, as reshaping simulations present challenges to be solved as the elevated computational cost of each simulation or the use of static data of traditional Finite Element Analysis (FEA), to name a few.

This report describes the distortions on aluminum forgings, the operations required to solve this problem and the challenges associated with deploying numerical simulations at workshop level and actually help operators. The document is structured as follows: Section 2 provides an overview of distortions as an open problem for the aeronautical industry. Then, the relation between residual stresses, distortions, and reshaping are explained in Section 3. Next, the different reshaping operations are described in detail in Section 4 and the challenges and perspectives to simulate reshaping are reviewed in Section 5. Finally, the obtained conclusions are presented.

Distortion in the aeronautical industry: an open problem

To the best of our knowledge, the problem of distortions in the aeronautical industry was first reported during World War II, when some difficulties were experienced with extrusions used for aircraft spars parts [START_REF] Forrest | Internal or Residual Stresses in Wrought Aluminium Alloys and their Structural Significance[END_REF]. Then, due to the technology transfer performed after the war, the use of large forgings for light metals alloys (i.e magnesium and aluminum) was potentiated with governmental initiatives as the Air Force Heavy Press Program in the United States, which allowed the construction of forgings up to 50k tons capacity1 [START_REF]ASME: The Wyman-Gordon 50[END_REF]. Thanks to this new range of tools, the unitization 2 of structural parts was feasible and acted as a driving force for the development of bigger aircraft with important manufacturing costs savings [START_REF] Ball | Inclusion of Forging Residual Stresses in Large Component Structural Design[END_REF][START_REF] Smith | The Production of Large Forgings in Aluminium Alloys[END_REF]. On the other hand, as the forgings grew in size and with them their cross sections and wall thickness, the range of residual stresses developed after quenching experienced an increment too, so as the research interest to predict the residual stresses and any future distortion after machining [START_REF] Jeanmart | Finite element calculation and measurement of thermal stresses in quenched plates of high-strength 7075 aluminium alloy[END_REF][START_REF] Robinson | 50th anniversary article: The origin and management of residual stress in heat-treatable aluminium alloys[END_REF]. At an industrial level, different attempts have been made to minimize distortion along the manufacturing chain of heat-treated aluminum alloys components. In North America, the MAI (Metals Affordability Initiative) program [START_REF] Ball | Inclusion of Forging Residual Stresses in Large Component Structural Design[END_REF][START_REF] Ma | Modeling of Residual Stress and Machining Distortion in Aerospace Components[END_REF] and in Europe, the COMPACT (A COncurrent approach to Manufacturing induced Part distortion in Aerospace ComponenTs) project [START_REF] Robinson | Residual stress in 7449 aluminium alloy forgings[END_REF][START_REF] Sim | Challenges of residual stress and part distortion in the civil airframe industry[END_REF][START_REF] Sim | Residual Stress Engineering in Manufacture of Aerospace Structural Parts[END_REF], focused on the understanding of residual stress and part distortion. However, both initiatives efforts were concentrated on the upstream manufacturing processes, such as quenching, stress relief, ageing, and machining, letting reshaping unattended and as a byproduct. Theoretically, if distortion is controlled or minimized, then reshaping is unnecessary. Unfortunately, a part with minimized distortions is only achieved in exceptional cases [START_REF] Mitze | Straightening heat-treated components[END_REF]. Therefore, we can state that distortion is an open problem for the industry during the last 80 years.

Residual Stresses, Distortion and Reshaping

In the distortion mitigation problem, there are three main actors: residual stresses, distortion, and reshaping. They can be considered as the source, the problem, and the current solution, respectively. At the same time, they are interconnected and interact between them, as represented in Figure 2. In the following section, a definition for each actor is provided.

Residual stresses refer to any stress distribution, which is present in a solid in the absence of an external load or thermal gradient [START_REF] Schijve | Residual Stress[END_REF]. In a material, the residual stress distribution is a direct consequence of non-homogeneous plastic defor- mations arising from transformations endured during the manufacturing process. Based on its characteristic length scale, residual stresses are classified as type I, II and III, corresponding to macro stresses, intergranular stresses, and stresses at an atomic scale, respectively [START_REF] Withers | Residual stress. Part 1 Measurement techniques[END_REF]. For distortion generation/mitigation, only residual stresses of type I are considered as this scale includes the stresses along the cross-section. This type of stresses can be induced by thermal or mechanical loading (e.g. welding and/or ultrasonic shoot penning).

Sim [START_REF] Sim | Challenges of residual stress and part distortion in the civil airframe industry[END_REF] defines distortion after machining as the shape deviation respect to the design geometry once the part is released from the fixture clamps. This problem is not caused by machining tolerance, over/under machining or dimensional inaccuracy. In a general way, distortion can be classified into two categories: in-plane and out-of-plane. Nevertheless, more complex patterns are present in forgings parts. Marin [START_REF] Marin | Calcul et optimisation des structures mécaniques[END_REF] documented the principal types of distortions present on a cruciform as shown in Figure 3, where each type of distortion is represented in an isolated way. However, in a real part, all five types of distortions are combined simultaneously with the additional challenge that every time a new part is machined, a different magnitude of distortion is obtained.

Reshaping can be defined as the post-machining operations required to mitigate distortion. The general idea is to induce distortion on the contrary direction of the actual warped geometry to counteract the geometrical defect and restore the design geometry within a given tolerance. To perform this task, the material must be taken to the plastic regimen with the aim to induce plastic strains in the part. A detailed description of the reshaping operations is provided in the next section. Now that the main actors are defined, it is easier to explain the whole picture of distortion. First, residual stresses appear after quenching. This heat treatment is the main source of residual stresses as the core of the part and its surface are subjected to different and strong thermal gradients. As a consequence, once the Fig. 3. Principal types of distortion present on an A320 cruciform after machining. Adapted from [START_REF] Marin | Calcul et optimisation des structures mécaniques[END_REF] part reaches the room temperature, there is a bi-axial compression stress distribution at the surface which is counterbalanced by a three-dimensional tension state in the core [START_REF] Robinson | Residual stress in 7449 aluminium alloy forgings[END_REF]. However, its internal distribution evolves each time a new operation is performed in the manufacturing chain. After quenching, the part continues with a stress relief process, followed by a second heat treatment called ageing, both processes aiming at minimizing residual stresses before machining. Then, during machining, the material removal breaks the equilibrium state of the current residual stresses, producing a redistribution of residual stresses. As a consequence, distortion appears. Next, each time reshaping is applied, new residual stresses are introduced to the part until the theoretical geometry is restored within the design's tolerance.

Post-machining distortion mitigation techniques

Once distortions are present after machining, it is mandatory to remove the geometrical deviations before proceeding to the assembly stage. In that way, the installation of the structural part into the corresponding sub-component is easier and the addition of unexpected stresses into the system is avoided. To perform this task, there are four reshaping operations available: bending straightening, torsion straightening, roller burnishing, and ultrasonic needle peening, as shown in Figure 4.

Bending straightening is an elastic-plastic operation. Depending on the part's geometry and dimensions, it can be used to remove globally a longitudinal distortion or locally the wing's waving (see Figure 3). After unloading, the residual stresses along the cross-section present the so-called three zero-crossing distribution, as depicted in Figure 5a. From a material science point of view, bending straightening and its non-homogeneous deformation produces a rearrangement of existing residual stresses [START_REF] Ellermann | Residual stress states as a result of bending and straightening processes of steels in different heat treatment conditions[END_REF], causes non-homogeneous material properties, anisotropy of strength and Bauschinger effect [START_REF] Schott | Consequences of bend straightening processes on residual stresses and strength of quenched and tempered steels[END_REF].

Similar to the previous process, torsion straightening is an elastic-plastic operation too. For simple cylindrical geometries, after unloading, the residual stresses present the so-called three zero-crossing distribution, as depicted on Figure 5b but instead of the compression-tension stresses generated during bending, shear stresses are produced. For more complex geometries as rectangular or I sections, the membrane analogy can be used [START_REF] Timoshenko | Strength of Materials: Part II Advanced Theory and Problems[END_REF] and it is reported that the plastic range is obtained first at the round fillets for I beams [START_REF] Christopherson | A Theoretical Investigation of Plastic Torsion in an I-Beam[END_REF][START_REF] Pi | Inelastic Torsion of Steel I-Beams[END_REF]. This operation is used to correct any twisting distortion and as its bending counterpart, suffers from the presence of springback during unloading which has to be taken into account.

Roller burnishing is a surface finishing method that induces a plastic deformation under cold working conditions by applying pressure through a cylindrical tool [START_REF] Murthy | Burnishing of metallic surfaces -a review[END_REF]. In the context of reshaping, it is used to modify the distribution of residual stresses by inducing a local compressive residual stress at the surface and for our specific case, two rollers are used, as shown in Figure 4c. This configuration produces a symmetric residual stress distribution as depicted in Figure 5c. By applying this method, surface irregularities are removed. As a consequence, a smooth and shiny surface with better fatigue strength is obtained [START_REF] Akkurt | Comparison of roller burnishing method with other hole surface finishing processes applied on aisi 304 austenitic stainless steel[END_REF]. The depth of the work-hardening layer depends on the process parameters and material, but in general, it is small and varies between 0.2-1 mm [START_REF] Coules | Threedimensional mapping of the residual stress field in a locally-rolled aluminium alloy specimen[END_REF][START_REF] Murthy | Burnishing of metallic surfaces -a review[END_REF].

Ultrasonic needle peening is a special case of shot peening. Here, instead of using the traditional spheres to hit the part and induce compressive residual [START_REF] Akkurt | Comparison of roller burnishing method with other hole surface finishing processes applied on aisi 304 austenitic stainless steel[END_REF][START_REF] Forrest | Internal or Residual Stresses in Wrought Aluminium Alloys and their Structural Significance[END_REF] stresses at the surface of the component, oscillating cylindrical bars in the ultrasonic range (over 18 kHz) are employed [START_REF] Yin | Overview of ultrasonic shot peening[END_REF]. This operation stretches plastically the area where is applied and, in order to counterbalance this effect, it produces an S-shape residual stress profile starting with a compressive layer followed by a tension zone underneath, as represented in Figure 5d [START_REF] Schijve | Residual Stress[END_REF]. The method is characterized to be manual, fast and easy to apply, as the tool can be easily handled by the operator. On the other hand, the area of influence of the residual stresses is smaller respect to the cross-section and normally, it is applied to a single side of the part, giving as a result an asymmetric residual stress distribution.

The corrective stage of reshaping is an iterative operation and its work-flow diagram is represented in Figure 6. Starting from the initial geometry, a measurement and control step is performed in order to determine the level of distortion. Based on the result, one of the four reshaping operations described before is selected. To perform one single plastic deformation step, the selected reshaping operation needs to be configured, it means, the technological parameters µ i (e.g position, imposed displacement, length of the application) need to be defined. This process is repeated until the part is repaired and fits into the design tolerance.

Reshaping simulation: challenges and perspectives

As stated before, nowadays, reshaping relies on the experience of a group of welltrained operators, the boilermakers. The goal to introduce numerical simulation in this process is to assist the operators and provide support during the decisionmaking step in order to reduce the number of operations needed to repair a part. However, from a technical point of view, reshaping simulation presents the following challenges:

the initial distortion after machining δ 0 is unknown and, with the available state-of-the-art simulation techniques, there is a mismatch between the numerical prediction and the real distortion.

distortions vary from one part to other (non-reproducibility), that is why reshaping cannot be standardized. simple 2D calculations are not able to represent properly neither the parts geometry nor the distortion phenomenon, and as a consequence, 3D simulations are required.

each time a new reshaping operation is performed, distortion changes and with it, the next reshaping operation and its technological parameters. to implement the numerical assistance in an industrial context, a real-time response is required.

The first challenge is related to an accurate prediction of residual stresses along the manufacturing chain. Theoretically, if the residual stress in a forged part is known before the machining step, it would be possible to determine the best location for the final geometry inside the blank, which minimizes any possible distortion. This approach works well for parts made from rolled plates [START_REF] Cerutti | Methodology for aluminium part machining quality improvement considering mechanical properties and process conditions[END_REF][START_REF] Chantzis | An industrial workflow to minimise part distortion for machining of large monolithic components in aerospace industry[END_REF] due to the residual stress field is easier to determine compared to forged parts [START_REF] Sim | Challenges of residual stress and part distortion in the civil airframe industry[END_REF]. For the latter, important deviations from the numerical prediction respect to the real distortion where reported by Cerutti [START_REF] Cerutti | Numerical modelling and mechanical analysis of the machining of large aeronautical parts: Machining quality improvement[END_REF], as shown in Figure 7. As a consequence, any mismatch between the numerical and real initial distortion can invalidate any further reshaping simulation step.

The second limitation is due to the deterministic nature of numerical simulations. Once the solution has converged, if the algorithm is executed again, it is expected to obtain the same result (e.g. δ 0 ) while in reality, every time a new forged part is machined, although all the process parameters were not modified, a different value for δ 0 is obtained. This pattern is observed between parts of the same class and among symmetric parts (e.g left and right cruciform). Therefore, reshaping simulation needs to provide multiple possible initial distortion scenarios.

The use of 3D simulations to represent properly the distorted geometry implies the generation of numerical models with large degrees of freedom. It is a common practice to take advantage of any symmetry plane in order to simulate only one portion of the total part. Maybe the symmetry trick can be exploited during the first reshaping step. However, as the problem evolves and distortion gets more localized, the way each reshaping operation is applied obliges to simulate the whole structure.

About the fourth challenge, for each reshaping step i=1,2,3,..., the goal is to define the parameters associated with the operation µ i that minimize the remaining distortion δ i . It means that, for each reshaping step, it is required to run a multi-parametric analysis. The problem's complexity grows as previous and subsequent operations influence each other. Therefore, any error generated during a defined reshaping step will be forward propagated and can invalidate the simulation results.

In order to implement this numerical assistance tool in an industrial environment, it is required to access the simulation outputs in real-time, as the goal is to evaluate beforehand different scenarios for the reshaping step i and only then, to launch the real operation with the optimum configuration defined in advance.

Under the described scenarios, the requirement of multiple reshaping simulations together with adaptive modeling and a real-time response suggests launching a multi-parametric analysis under a Model Order Reduction (MOR) framework. The idea for this approach is to produce a low-dimensional but accurate model in a cost-effective way which allows exploring beforehand the input space to construct a solution dictionary, known as computational vademecum [START_REF] Chinesta | PGD-Based Computational Vademecum for Efficient Design, Optimization and Control[END_REF].

Conclusions

Distortion of large aluminum forgings is an open problem for the aeronautical industry. Continuous research efforts are done in order to predict the residual stresses that generate this problem. However, it is important from an operational standpoint to study how to handle distortion once it has appeared. In this document, the required post-machining operations to mitigate distortion has been reviewed together with their associated challenges to simulate them. At the moment, reshaping is solved based on the empirical experience acquired for a group of well-trained operators. Now, the goal is to find a way to use advanced numerical simulations tools (e.g adaptive modeling, model order reduction) to study reshaping in a virtual environment and evaluate possible scenarios before launching the real operation. The long-term goal is to optimize the process in order to minimize the total number of reshaping operations and avoid introducing unnecessary residual stresses to the structural part.

Fig. 1 .

 1 Fig. 1. Example of aluminum forgings subjected to distortion: cruciform of the A320 wingbox (depicted in red)

Fig. 2 .

 2 Fig. 2. Interaction between Residual Stresses, Distortion and Reshaping.

Fig. 4 .

 4 Fig. 4. Reshaping operations used at Airbus: (a) Bending, (b) Torsion, (c) Roller burnishing and (d) Ultrasonic needle peening.

Fig. 5 .

 5 Fig. 5. Residual stress process signature after unloading: (a) Bending, (b) Torsion, (c) Roller burnishing and (d) Ultrasonic needle peening. Adapted from[START_REF] Akkurt | Comparison of roller burnishing method with other hole surface finishing processes applied on aisi 304 austenitic stainless steel[END_REF][START_REF] Forrest | Internal or Residual Stresses in Wrought Aluminium Alloys and their Structural Significance[END_REF]

Fig. 6 .

 6 Fig. 6. Reshaping simulation workflow diagram

Fig. 7 .

 7 Fig. 7. Example of distortion in a cruciform. (a) Geometry before and after machining. (b) Simulation vs measurements results. Adapted from [5].

Nowadays, the biggest forging has a capacity of 80k ton and is located in China[START_REF] Du | Recent Progess of Manufacturing Technologies on C&W Superalloys in China[END_REF].

Acknowledgements. This project is part of the Marie Sk lodowska-Curie ITN-ETN AdMoRe funded by the European Union Horizon 2020 research and innovation program with grant number 675919.


Ramiro Mena

José V Aguado

Stéphane Guinard

Antonio Huerta

Post-machining distortion mitigation via bending straightening, a multiparametric ROM study

Large and thick aeronautical structural parts frequently exhibit significant distortions after machining because of the residual stresses generated along the manufacturing process. A post-machining shape correction phase, called reshaping, has to be applied to each part before proceeding to the assembly stage. At the moment this process relies entirely on the know-how and experience of a specialized operator. As a consequence, the process is expensive and time-consuming. Therefore, there is a need at an industrial level to solve this problem with the support of numerical simulation tools. However, to reach this target, it is required to explore beforehand the input space to construct a solution dictionary in order to provide real-time response. Both requirements oblige to launch a multi-parametric analysis under a Model Order Reduction (MOR) framework. A complete reshaping process consists of a sequence of elementary operations, among which bending straightening is selected in its three-point bending configuration and it is simulated by applying the Sparse Subspace Learning method (SSL). By using the MOR approach, new insights for the distortion mitigation via plastic bending can be explored that will help to understand the open problem of reshaping.

INTRODUCTION

Large and thick aeronautical structural parts are made of aluminum forgings in order to produce a lightweight construction [1]. However, they frequently exhibit significant distortions after machining because of the residual stresses generated during the sequential manufacturing steps, especially after the heat treatment of quenching [START_REF] Ma | Modeling of Residual Stress and Machining Distortion in Aerospace Components[END_REF][START_REF] Robinson | [END_REF]. Such distortions are variable in nature, but the industrial process incorporates a robust reshaping stage: skilled boilermakers operate series of mechanical loadings on the non-conforming part and re-establish its nominal geometry while obeying guidelines imposed by stress offices -which guarantees the produced part has still mechanical characteristics compatible with its operational role -. Nonetheless, numerical simulations may offer here huge margins for improvements: operational assistance to boilermakers and related gains in time-cycle (optimized reshaping sequence, optimized parameters for each elementary reshaping operation), improved knowledge of actual part and related gains in safety (monitoring of actual plastic strains) [START_REF] Sim | Challenges of residual stress and part distortion in the civil airframe industry[END_REF][START_REF] Mena | Influence of the residual stresses in reshaping operations of large aeronautical parts[END_REF]. An example of the structural parts subjected to reshaping are the cruciform beams that form the wing-box in an aircraft.

Reshaping can be considered as the sum of four local operations where the part is plastically deformed. These operations are: bending straightening, torsion straightening, roller burnishing, and ultrasonic needle peening. The main idea in this process is to produce a local deformation in the opposite direction of the deformed shape in order to counteract the mismatch between the real and the theoretical geometry.

This corrective stage is an iterative operation and its work-flow diagram is represented in Figure 1. Starting from the initial geometry, a measurement and control step is performed in order to determine the level of distortion. Based on the result, one of the four reshaping operations described before is selected. To perform one single plastic deformation step, the selected reshaping operation needs to be configured, it means, the technological parameters µ i (e.g position, imposed displacement, length of the application) need to be defined. This process is repeated until the part is repaired and fits into the design tolerance. As stated before, nowadays, reshaping relies on the experience of a well-trained operator. The goal to introduce numerical simulation in this process is to assist the operators and provide support during the decision-making step in order to reduce the number of operations needed to repair a part. However, from the technical point of view, reshaping simulation presents three main challenges:

• the initial distortion after machining is unknown and, with the available state-of-the-art simulation techniques, there is a mismatch between the numerical and the real distortion [START_REF] Cerutti | Numerical modelling and mechanical analysis of the machining of large aeronautical parts : Machining quality improvement[END_REF].

• distortions vary from one part to the other (non reproducibility), that is why reshaping cannot be standardized.

As a consequence, reshaping simulation needs to provide multiple possible initial distortion scenarios. • each time a new reshaping operation is performed, distortion changes and with it, the technological parameters associated with the next reshaping operation. Hence, reshaping simulation must be able to handle multiple parameters in real-time.

Under the described scenarios, the requirement of multiple reshaping simulations together with a real-time response obliges to launch a multiparametric analysis under a MOR framework [START_REF] Chinesta | [END_REF]. For the different MOR available techniques, the Sparse Subspace Learning (SSL) is the selected method to study reshaping thanks its non-intrusiveness and the capability to handle the nonlinearities present during the process (e.g geometrical non-linearity and plasticity). For a detailed description of the method, the interesed reader can be referred to [8,9]. This document is organized as follows: first, the SSL formulation applied to reshaping is explained. Then, the problem setup used for a T shaped beam geometry is described. Next, the results of the multi-parametric ROM study are analysed and finally, the obtained conclusions are addressed.

RESHAPING UNDER THE SSL FORMULATION

Sparse Subspace Learning (SSL) is a non-intrusive Model Order Reduction technique, suitable for multiparametric studies and nonlinear analysis. The method exploits the sparsity of the solution with respect to the parameters in order to build a parametric solution from the computation of only few wisely chosen snapshots [8]. From the four available reshaping techniques, the bending straightening operation is selected in its three-point bending configuration to be studied with the SSL method, due to this operation is the most used to repair the structural parts in the industry [10].

For bending straightening, we are interested in determining the final shape of a part after performing a cycle of loading and unloading. Therefore, the displacement field u(x) for a given set of parameters µ is our objective function and, by using the separated variable representation under the SSL formulation, it can be written as:

u (x; µ) = k i=0 α i (x) ψ i (µ) (1)
where α(x) and ψ(µ) are defined as the surpluses and parametric functions, respectively. The parametric functions ψ(µ) act as the shape function. They are known, define the hierarchical basis and correspond to the Lagrangian polynomials, while the surpluses functions α(x) need to be determined. They are calculated as a linear combination of sampling or measurements evaluated at the collocation points P n = {µ 0 , µ 1 , ...., µ n }.

For a unidimensional parametric problem, the collocation points are the extreme of the n-th Chebyshev polynomials of the first kind, T n (µ) and are referred as Gauss-Chebyshev-Lobatto points (GCL) too.

In order to handle a multidimensional parametric problem, first, each parameter is treated as the unidimensional case described before and then, the Smolyak rule is introduced with the aim of constructing the sparse grid to be used during sampling. To put into perspective the advantage of the sparse sampling respect to the full sampling based on the tensor product, for a 2 parameters problem and a selected hierarchy level k = 5, the full sampling requires 33 2 = 1089 points, while the sparse sampling only needs 177 points, providing a reduction of 83.7% in terms of computational cost.

Once the sampling points are defined, the simulations are launched in a commercial FEM software (e.g ABAQUS, PamCrash). The displacement field u (x; µ i ) is known at each collocation point and for our case, it is considered as the reference solution. With this available information, the corresponding surplus functions are determined as the difference between the exact and the approximated solution at the hierarchy level (k -1), as expressed in the following expression:

α i (x) = u (x; µ i ) -u k-1 (x; µ i ) (2)
Coming back to the reshaping problem, for the three-point bending straightening operation, the set of parameters µ is composed by the initial distortion δ 0 , the position of the bottom roller r b and the imposed vertical displacement y s . For our problem, each parameter is considered as a dimension. µ = (δ 0 , r b , y s )

When the solution u(x; µ) is reconstructed for all the parametric domain, the remaining distortion δ can be found as a post-processing operation and it is defined as:

δ (µ) = max x u (x; µ) -min x u (x; µ) (4)
However, for all the possible reshaping configurations, we are mainly interested to know the optimum stroke y opt s that, for a fixed set of parameters (δ 0 , r b ) will minimize the remaining distortion δ in the structural part.

y opt s = arg min y s δ (δ 0 , r b , y s ) (5)

PROBLEM SETUP

To study the three-point bending operation, a T shaped beam is defined, as shown in Figure 2. This geometry corresponds to a simplified version of a cruciform part, which exhibits critical distortions when machined and is a primary structure of main importance for aircraft manufacturers. The parameters to be explored are available in Table 1 with is minimum and maximum values µ min and µ max , respectively. The selected material is an aluminum AA7010 previously characterized by Airbus Group Innovations. The elastic material properties are summarized in Table 2. To run the reshaping simulations, a Chaboche's isotropic/nonlinear kinematic hardening model was used. However, as the material parameters are under a nondisclosure agreement, they are not included in this document.

Finally, a Python script in combination with ABAQUS was used to run the 123 required simulations to explore the parametric input space under the SSL method for a hierarchy level k=3.

RESULTS AND DISCUSSION

Before proceeding to explain the effect of each parameter in the bending straightening operation, it is required to introduce the evolution of the remaining distortion δ along the process. An example of the different zones present during reshaping is represented in Figure 3. In the beginning, if a small value of y s is applied, the warped geometry will remain inside of the elastic region and there is no modification of δ 0 . We define this region as zone A and it is characterized by the relation δ ≈ δ 0 . A reduction of the initial distortion is only performed when plasticity takes place along the part. We call zone B at the region where δ < δ 0 . Inside this zone, the optimum configuration (δ = δ min ) is located in a narrow space. To find properly this configuration is one of the challenges of reshaping. On the other hand, an excessive value of y s can be counter-productive and invert the shape of the initial distortion from concave to convex (or vice-versa), giving as a result the operation zone C, where δ > δ min , and which should be avoided. The reshaping diagram described before is the result of a particularized solution where δ 0 and r b were previously defined and y s was set free. In order to generalize the problem, the next step is to set free a pair of parameters (e.g δ 0 and y s ), giving, as a result, a response surface. With the aim to explore the obtained multi-parametric solution, instead of using a 3D representation, a 2D plot together with level curves are used to describe if the remaining distortion δ will decrease or increase for a fixed bottom roller position r b , as can be appreciated in Figure 4.

Effect of the initial distortion δ 0 : In a previous work, it has been proved that reshaping can be studied by using as a main input δ 0 under a Residual Stress Free hypothesis (RSF) [START_REF] Mena | Influence of the residual stresses in reshaping operations of large aeronautical parts[END_REF]. This approach was used in order to generate numerically the warped geometry. By including δ 0 as a parameter, the problem of not knowing its magnitude beforehand is overcome. In that way, the reshaping problem is solved numerically for all the interval of initial distortions and once the part is machined, the operator only needs to access to the solution dictionary and select a distorted geometry closer to the experimental result.

Effect of the bottom roller position r b : Once the initial distortion δ 0 is known, the next step is to select the reshaping operation. For the three-point bending operation, the internal roller (top) is fixed and located at the middle of the beam while the external rollers (bottom) need to be placed. The position of the bottom roller is defined by the parameter r b and it defines the slope of the response surface that relates the remaining distortion δ with δ 0 and y s .

A value of r b closer to the top roller produces a difficult reshaping operation to control, as it is expected to have a quick change of the remaining distortion δ during bending. On the other hand, a value of r b away from the top roller generates the opposite behavior with a progressive variation of δ. In the performed parametric study, the selected value of r b particularizes a page of the solution dictionary. Effect of the imposed vertical displacement y s : For a fixed initial distortion and a selected bottom roller position, the three-point bending operation is completed when a vertical displacement y s is imposed at the internal roller. The magnitude of y s needs to be defined in order to minimize the remaining distortion δ. The response surface obtained for a fixed r b shows that reshaping is more difficult to be performed for smaller values of δ 0 as the required stroke y s only to induce plasticity in the cross section is bigger than δ 0 . On the other hand, for a bigger magnitude of δ 0 , it is expected to have more margin to find the stroke that minimizes δ.

CONCLUSIONS

The industry is willing to have new tools to solve problems that nowadays are treated in an empirical way. Reshaping of large thick structural parts is an example of how MOR can be applied to an open problem, especially the SSL, allows to perform a multi-parametric study to test different what-if scenarios in a cost-effective way and, thanks to its non-intrusiveness, acts as a bridge to join MOR with the available third-party FEM solutions.

For the case of bending straightening, the initial distortion δ 0 is the starting point for reshaping and it determines the setup for the technological parameters r b and y s . Additionally, the problem of not knowing beforehand the value of δ 0 is overcome by including it as a parameter. Then, thanks to the SSL, this reshaping operation can be studied in a virtual environment before launching the real process.

Future work will focus on how to chain multiple reshaping operations and to include material properties as extra parameters. As expected, the computational cost for the numerical analysis will increase. Again, here it resides the importance of studying reshaping with MOR techniques.

FIGURE 1 .

 1 FIGURE 1. Reshaping work-flow.

FIGURE 2 .

 2 FIGURE 2. T shaped beam. a) Geometry definition and b) Parameters representation for three-point bending straightening. TABLE 1. Parameters values for the three-point bending operation Parameter ID Description µ min µ max 1 Initial distortion δ 0 (mm) -15 0 2 Bottom roller position r b (mm) 100 485 3 Imposed vertical displacement y s (mm) 5 30

FIGURE 3 .

 3 FIGURE 3. Zones present during reshaping as a function of the stroke y s .

TABLE 2 .

 2 Elastic properties for Aluminum AA7010

	73000	0,3	390

Young's Modulus E (MPa) Poisson's ratio ν (-) Yield stress σ y (MPa)

ACKNOWLEDGMENTS

This project is part of the Marie Skłodowska-Curie ITN-ETN AdMoRe funded by the European Union Horizon 2020 research and innovation program with grant number 675919.


Ramiro Mena

José V Aguado

Stéphane Guinard

Antonio Huerta

Reshaping of Large Aeronautical Structural Parts: A Simplified Simulation Approach

Large aeronautical structural parts present important distortions after machining. This problem is caused by the presence of residual stresses, which are developed during previous manufacturing steps (quenching). Before being put into service, the nominal geometry is restored by means of mechanical methods. This operation is called reshaping and exclusively depends on the skills of a well-trained and experienced operator. Moreover, this procedure is time consuming and nowadays, it is only based on a trial and error approach. Therefore, there is a need at industrial level to solve this problem with the support of numerical simulation tools. By using a simplification hypothesis, it was found that the springback phenomenon behaves linearly and it allows developing a strategy to implement reshaping at an industrial level.

INTRODUCTION

In the aircraft industry, the technical requirement of minimum weight with maximum strength generates parts with complex geometries. As a result, large monolithic structures are made of aluminum forgings. These components present residual stresses as a consequence of non-uniform plastic strains generated during previous manufacturing steps, especially after quenching 1,2 . When parts are machined, in order to obtain the final geometry, almost 90% of material is removed under high feeds and high cutting speeds. This material removal brakes the previous stable configuration and a redistribution of residual stresses take place 3 . As a result, important distortions that are out of the geometrical tolerances arise. Before being put into service, distortions are removed manually by mechanical methods. This operation is called reshaping and exclusively depends on the skills of a well-trained and experienced operator. Moreover, this procedure is time consuming and nowadays, it is only based on a trial and error approach. Therefore, there is a need at industrial level to solve this problem with the support of numerical simulation tools [START_REF] Sim | Challenges of Residual Stress and Part Distortion in the civil airframe industry[END_REF] .

Reshaping can be interpreted as an optimization problem, where the goal is to introduce the minimum number of extra reshaping steps in order to fit the geometry inside a given tolerance. The idea behind this post-machining operation is to produce localized plastic strains in determined part's zones to minimize remaining distortion. One common reshaping operation consists in loading the distorted part backward (bending in the opposite direction to observed distortions): during unloading, as the external load is removed, the elastic component of strain is released and the springback phenomenon is present. As a consequence, an extra plastic strain is required to compensate this issue. Although it was originally applied to straightening of seam pipes, Katoh et al., proposed to approximate the real springback behavior by a linear hypothesis in a force-displacement diagram [START_REF] Katoh | [END_REF] . This idea is translated to our problem and tested for two different geometries: a rectangular beam and a T shaped beam. The first case is solved experimentally to validate the range of applicability of Katoh's hypothesis, whereas the second case is solved numerically. This last model corresponds to a simplified version of a cruciform part, which exhibits critical distortions when machined and is a primary structure of main importance for aircraft manufacturers. This document is organized as follows: first, the idea of the force displacement diagram applied to reshaping is explained. Then, the problem setup used for the two beam geometries is described. Next, the results of both experiments are analyzed and finally, the obtained conclusions are addressed.

FORCE DISPLACEMENT DIAGRAM

Given an initial distortion, a warped part is straightened by using mechanical operations, normally a three or four point bending. The first case is selected for localized distortions such as kinks whereas the last is used when distortion is progressively distributed along the part. For the rest of the study, the four point bending is chosen and it is represented in Fig. 1 (a), where L defines the effective span length, δ r corresponds to the reference distortion and a is the distance between roller supports and the applied force F.

For a defined technological setup (i.e roller positioning), there is an imposed given displacement, called stroke Y s , capable to minimize distortion, where the springback δ s is properly compensated. In order to find this configuration, a traditional FEM algorithm obliges to run a new simulation each time Y s is modified. A more suitable approach, proposed by Katoh et al., used the force F -displacement δ diagram to compensate in real-time the springback behavior by assuming that the gradients during loading and unloading are parallel [START_REF] Katoh | [END_REF] . It means that the equivalent Young's modulus E 0 , calculated during loading, can be used for unloading Ê 0 , E 0 ≈ Ê 0 . Therefore, as δ r is an initial data and corresponds to the distortion to be removed, the loading stage during reshaping will take place until the resultant F-δ curve intersects the predicted unloading path.

In the F-δ diagram (Fig. 1 (b)), during loading, the slope of the curve is defined by E 0 , and follows the relation:

δ 0 E F = (1)
For the selected reshaping operation, E 0 can be found theoretically by using the beam expression:

= 24 / [ 3 -4 ] (

) 2
where E and I are the Young's modulus and inertia of the cross section, respectively. Another option is, given two load and displacement points before yielding, named as (F 1 ,δ 1 ) and (F 2 , δ 2 ), respectively, E 0 can be determined experimentally as:

( ) ( ) 1 2 1 2 0 δ δ - - = F F E (3)
this last approach has the advantage that no previous material modeling information is needed. Additionally, when it comes to a part with residual stresses, it allows overcoming the uncertainties associated with simulating each manufacturing step, simplifying the problem and adapting the response to each specific part.

PROBLEM SETUP

In order to study reshaping, two models are proposed. First, starting from a flat configuration, a rectangular beam with chamfers at the top surface is plastically deformed by imposing a vertical displacement of ±9,5mm during three cycles. The sample's geometry and the test setup are represented at Fig. 2 (a) and (b), respectively. For this test, displacements are recorded at two locations: at the beam's mid-span and where the force is applied. Additionally, the force data acquisition is done by the hydraulic press' head and strain data at the top and bottom faces are recorded by two strain gauges, located at the intersection of both longitudinal and transversal symmetry planes of the beam. The goal of this test is to validate the applicability of Katoh's hypothesis for a simple geometry subjected to a deformation in the range of 4%.

The second model corresponds to a T shaped beam. This case is based on the results previously obtained in reference [START_REF] Mena | Influence of the residual stresses in reshaping operations of large aeronautical parts[END_REF] , where distortion was generated numerically and two configurations were selected, named as configuration A and B, respectively (see Fig. 3 (c)). Then, reshaping was done by incrementing sequentially the imposed vertical stroke Ys until distortion was minimized. The results were calculated by using a nonlinear material behavior. Now, they are compared by assuming an unloading linear path.

For both experiments, the selected material is aluminum AA7010, which was previously characterized by Airbus. Its elastic properties are summarized in Table 1. To run the reshaping simulations, a Chaboche's isotropic/nonlinear kinematic hardening model was used. However, as the material parameters are under a nondisclosure agreement, they are not included in this paper.

RESULTS AND DISCUSSION

For the first experiment, the force-displacement diagram is obtained at two locations: the beam's mid-span and the section where the imposed displacement is applied. Both diagrams are represented in Fig. 4 (a) and (b), respectively. The equivalent Young's modulus E 0 is calculated by using equation (3), then, a linear unloading path is represented at the final position of each cycle. It can be seen that the experimental unloading path presents a small discontinuity at the beginning and then, it follows a nonlinear trajectory. However, this trajectory is convex and for reshaping, we are interested in two points: where to stop the loading curve and what is the final distortion in the part. Therefore, if the linear approximation of the unloading path is used, both points can be tracked in a straightforward way for a given tolerance.

An explanation to Katoh's hypothesis can be found in literature [START_REF] Yu | Plastic Bending: Theory and applications[END_REF] , where Yu et al., state "When a complete unloading is performed after applying a bending moment M to a beam, if unloading does not cause reverse yielding, then effectively the unloading process is equivalent to the elastic effect caused by applying -M to the beam". Here, the key is to not produce reverse yielding. This phenomenon is present in the automotive industry during metal sheet forming. However, the structural elements that are subjected to reshaping in the aircraft industry can be considered as thick walled, and the levels of strain under reshaping cannot produce reverse yielding. Therefore, the unloading path can be approximated to a linear behavior. A similar analogy can be done respect to the residual stresses generated during machining, which are required to be included in a distortion simulation for a thin walled component 8 or that can be neglected in a thick walled part [START_REF] Cerutti | Numerical modelling and mechanical analysis of the machining of large aeronautical pats: Machining quality improvements[END_REF] .

If the lineal hypothesis is valid during unloading, then superposition can be used and the residual stresses σ residual can be found by adding the stresses generated during loading σ load to the unloading stresses σ unload :

unload load residual σ σ σ + = (4)
By reproducing numerically the experiment, for the values obtained during the 1 st cycle, σ load and σ residual are known, which correspond to the stress states at the beginning and at the end of the unloading step, respectively. Then, σ unload is obtained by using expression (4). The previous stress states are shown in Fig. 5 and they are represented along the cross section height (y-axis) at the intersection of the longitudinal and transversal symmetry planes. Again, it can be seen how σ unload presents a nonlinear distribution (Fig. 5 (b)). When fitted linearly, an approximated residual stress distribution is obtained (Fig. 5 (c) in dashed lines), which provides a close approximation of the result obtained numerically. After each cycle, the residual stress (RS) follows a triple zerocross pattern [START_REF] Ellerman | [END_REF] , its order of magnitude is around ± 200 MPa after unloading and bend straightening causes RS similar to the values generated after quenching. For the second experiment, corresponding to the T shaped beam, the numerical results were obtained in a previous work with a nonlinear material model [START_REF] Mena | Influence of the residual stresses in reshaping operations of large aeronautical parts[END_REF] . A similar procedure as the one previously described is used to calculate E 0 . Again, a force-displacement diagram is generated for reshaping configurations A and B, which are represented in Fig. 6 (a) and (b), respectively. As it can be seen, the linear unloading path approximation is valid for both reshaping configurations, although, configuration B presents a more pronounced nonlinear unloading path, probably caused by the higher initial distortion. From the multiple available stroke scenarios, we are only interested in the parameter that minimizes distortion, which corresponds to the red curves depicted in both figures. For the optimum curves, if the linear unloading path is used, the difference in terms of maximum load during reshaping is less than 3%. The corresponding values for configuration A and B are summarized in Table 2. Additionally, if both configurations are considered as the extreme values for the interval of machining offset Mo = [5, 10]

CONCLUSIONS

Reshaping simulations can be considered computationally expensive, as 3D solid elements are used and multiple simulation steps are required before obtaining an initial warped geometry. Additionally, an adaptive response is desired to modify a simulation if some reshaping conditions are changed at the workshop level. Therefore, some simplifications hypotheses are needed in order to implement the numerical results in an industrial environment. Following this line, the idea of considering a distorted part free of residual stresses for simulating reshaping was able to describe the part's distortion evolution compared to the case with residual stresses [START_REF] Mena | Influence of the residual stresses in reshaping operations of large aeronautical parts[END_REF] , however, its resolution remains subjected to the progressively variation of the imposed vertical stroke Ys until distortion δ is minimized. Now, by using a linear approximation to the unloading path, only one loading simulation would be required and multiple unloading scenarios can be easily explored. Finally, the force and displacement measurement are parameters that can be implemented in an industrial environment.

FIGURE 1 .

 1 FIGURE 1. Reshaping operation. (a) Initial warped geometry, (b) Force vs displacement diagram.

FIGURE 2 .

 2 FIGURE 2. Rectangular beam experiment. (a) Geometry definition, (b) Laboratory setup.

FIGURE 3 .TABLE 1 .Young

 31 FIGURE 3. T shaped beam model. (a) Geometry definition, (b) Initial distortion as a function of machining offset M O , (c) Distorted geometry cases and reshaping setup. TABLE 1. Elastic properties for Aluminum AA7010 Young's Modulus E (MPa) Poisson's ratio ν (-) Yield stress σ y (MPa) 73000 0,3 390

FIGURE 4 .

 4 FIGURE 4. Force vs displacement diagram for the rectangular beam. Experimental curves (solid lines) vs linear unloading path (dashed lines). (a) Mid-spam L=100 and (b) Imposed stroke Ys=±9,5mm.

FIGURE 5 .FIGURE 6 .

 56 FIGURE 5. Longitudinal stress σ 11 for the rectangular beam. (a) Loading, (b) Unloading and (c) Residual

TABLE 2 .

 2 Starting point for the unloading curve at the optimum stroke (named as Opt in Fig.6)

	Configuration	Numerical max. load (KN)	Simplified max. load (KN)	Error (%)
	A	-25.99	-26.7	2.732
	B	85.5	86.52	1.193

ACKNOWLEDGMENTS

This project has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 675919.


J B Fl6r

J M Burgers

H J S Fernando

G J F Van Heijst

The evolution of an isolated turbulent region in a two-layer fluid

A turbulent region is geperated by horizontal pulsed injection at the interface of a two-layer fluid. Flow visualizatiod; studies reveal the existence of three stages in the evolution of the vertical size of this region: growth, maximum height, and collapse. A scaling analysis for the height of the turbulent region is presented, which appears to be in good agreement with the measurements. Comparable results were obtained by Fernando, van Heijst, and Fonseka (submitted to J. Fluid Mech.) for similar experiments in a linearly stratified fluid. Thorpe-scale measurements of the turbulent region reveal that the ratio of the rms displacement L1 and the maximum displacement Ltrnax. remain constant with time. The eventual formation process of the dipolar vortices after the collapse and the influence of interfacial wave motions on these dipolar vortices are discussed.

I. INTRODUCTION

In many geophysical flows stable stratification plays a major role in the evolution of turbulent flows. In such flows, vertical kinetic energy of fluid particles is converted into potential energy during their overturning motions, which leads to the suppression of vertical spreading of the turbulent region. After having reached its maximum ver tical size, the mixed region starts to physically collapse under gravity until the motion is degenerated into a hori zontal eddy pattern.

Local mixing events in stratifi ed flows can occur by the breakdown of Kelvin-Helmholtz billows; some examples are the oceanic thermocline observations of Woods 1 and Gregg, 2 and the atmospheric observations of Browning. 3 Other processes responsible for mixing events in stratified fluids are the breaking of internal waves, critical-layer ab sorption, and double diffusion. [START_REF] Marmorino | Observations of small-scale mixing processes in the seasonal thermocline. Part II: Wave breaking[END_REF] In the laboratory, the evolution of an isolated turbu lent region in a linearly density-stratifi ed fluid has been investigated by Fernando, van Heijst, and Fonseka, 5 here after referred to as FvHF. In that study, an isolated tur bulent region was generated by horizontally injecting a small volume of fluid during a short period of time, and the evolution of the emerging turbulent region was investi gated. In this evolution three different stages have been identified: the turbulent region grows vertically (first stage) until it reaches its maximum size (the second stage), and subsequently collapses under gravity (third stage). The flow eventually is transformed into a dipolar vortex structure, as was described by van Heijst and Fl6r6 and Voropayev et al. [START_REF] Voropayev | Horizontal jets and vortex dipoles in a stratified fluid[END_REF] For each stage a scaling analysis was given for the evolution of the vertical patch size.

The purpose of the present investigation was to extend the work of FvHF to a two-layer fluid; i.e., the pulsed injection of fluid takes place at the interface of a two-layer fluid. Although the details are different, the qualitative fea tures of the evolution in linearly stratifi ed and two-layer fluids appear to be the same. In the present study, flow visualization and measurements of the vertical density dis tributions are tested to study the evolution of the injected turbulent region. The experimental arrangement for these measurements is described in Sec. II. Qualitative observa tions are described in Sec. III, and the experimental results as well as scaling analyses of the growth and the subse quent collapse of the turbulent region are presented in Sec.

IV. Additional measurements on the dipolar vortex struc ture that is formed at the later stages are discussed in Sec.

V. A brief summary of the results and conclusions are given in Sec. VI.

II. EXPERIMENTAL SETUP

The experiments reported in this paper were conducted in two different Perspex tanks: one narrow tank of dimen sions 183X30X40 (em) and one square tank of dimen sions 90X90X40 (em), respectively. For the study of the growth and collapse of the turbulent patch, the narrow tank was used, while the square tank was used to study the dipole formation. In order to assure that side wall effects in the small tank did not bias the collapse significantly, the interface was siphoned off using small perforated horizon tal tubes. By this method the interface thickness oh was reduced to a thickness of approximately 0.6-2.5 em, which was always smaller than the size of the turbulent patch. The total depth of the water was varied from 20 to 31 em, and the position of the interface, measured from the bot tom of the tank, was varied between 12 and 14 em. The buoyancy jump across the interface, represented by

8b=2g p z -Pt ' P1+P2 (1)
with p 1 and p 2 the densities of the upper and lower layer, respectively, and g the gravitational acceleration, was in the range from 0.37 to 0.77 m/s 2 • In order to measure the density profile a conductivity probe was traversed verti cally through the fluid, from which the interfacial position was determined. Then the nozzle was placed at this level and a small drop of the dyed injection fluid-with density ( p 1 + p 2) /2-was released at the level of the nozzle. Due to the strong buoyancy force at the interface the drop col lapses instantly to a thin spot of less than a millimeter thickness and marked the midlevel of the interface quite accurately. The fluid was injected horizontally at this level through a nozzle of diameter 2 mm by a computer controlled injection mechanism. The density of the injected fluid was adjusted to the intermediate density, (p1 + p 2 )/2. The injection velocity in these experiments was varied from 1.8 to 11.7 m/s, and the injection volume V0 ranged from 1 to 13 ml. The time of injection 8t was varied between 0.01-0.2 s, which was sufficiently short to create an iso lated turbulent patch of fluid. The Reynolds number Re= U0d!v at the injection, where U0 is the injection ve locity, d is the nozzle diameter, and v is the kinematic viscosity, ranged from 3500 to 22 000, thus ensuring that the jet was turbulent. The resulting turbulent patch was visualized by coloring the injected fluid with fluorescein and by illuminating it using a vertical sheet of light. The extent of dye diffusion was considered as the extent to which the turbulence diffuses; that is, the turbulent Prandtl number is considered to be of the order unity (Refs. 8 and 9). A horizontal light sheet was used for visualizing the dipolar vortex formation in the square tank. The evolution of the flow was recorded by a video camera and by a frame by-frame examination of these recordings the time evolu tion of the vertical patch size was determined. The moment at which the fi rst spot of colored fl uid is ejected by the nozzle, was taken as time t=O. The vertical density distribution was measured by us ing a fast traversing microscale conductivity probe that was installed above the square tank. This traversing mech anism could be translated horizontally to track and capture the moving patch during the course of an experiment. The plunging speed of the probe had a maximum value of 30 cm/s, which is fast compared with the observed maximum translation speed of the patch (about 5 cm/s). These in stantaneous vertical density profi les were used to calculate Thorpe displacements, i.e., the size of the vertical displace ments of fluid parcels relative to a stable profile, which was obtained by rearranging the observed density profile. The method employed is similar to that used by Dillon. 10 The time and position at which the probes were shot, were recorded by video.

Ill. QUALITATIVE OBSERVATIONS

The evolution of the turbulent patch after the injection is presented in Fig. 2. During and shortly after the injec tion, an approximately cone-shaped jet forms and evolves into an isolated three-dimensional turbulent patch and moves in the direction of injection. During the initial stage, the vertical inertia forces of the eddies are large compared to the buoyancy forces, and entrainment of homogeneous fluid from above and below the interface causes the turbu lent region to grow continuously [see Fig. 2(b)]. Mixing of jet fluid and ambient fluid implies that the kinetic energy is converted into potential energy. When the vertical inertia forces become of the same order as the buoyancy forces, the turbulent patch reaches its maximum vertical size and then starts to collapse under gravity, thus reducing its ver tical size [for example, see Fig. 2(c)]. During the collapse, the mixed fluid intrudes horizontally into the ambient fluid. In some cases one could observe the formation of two separated intrusion noses at the front of the intrusion: one above and one below the interface [see, e.g., Figs. 2(c) and2(d)]. The occurrence of the double intrusion can be at tributed to incomplete mixing within the patch, due to the finite thickness of the initial interface, 11 as well as interfa cial friction between the two intrusions. The collapse con tinues until the patch attains a thickness comparable to that of the interface. Interfacial waves, generated during the collapse, could not be distinguished from the waves that were forced by the initial impulse of the jet. In most of the experiments it is observed that these waves move away from the turbulent region along the upper and lower , boundaries of the intrusion. These traveling waves are dis sipated rapidly, in particular, after having been reflected at the tank walls. From the video recordings, however, it ap pears that in some cases a soliton-like disturbance moves at the front of the patch, giving rise to a vertical bulge near the nose of the intrusion. This bulge vanishes slowly com- pared to the collapse rate of the patch. In a few experi ments with high injection velocities, a soliton-like feature was observed to detach from the collapsing region and travel up and down the tank a few times (see Fig. 3). Similar wave motions have been observed by Gilreath and Brandt 12 in the wake of an obstacle towed through a den sity interface. At the earlier stages of the collapse, the in trusion front is contorted and is observed to be distorted by flattened large eddies, which persist for longer times than the small-scale turbulence. At later stages, the front changes into one single, sharp intrusion. Mter the collapse the motion is mainly horizontal. A dipolar vortex structure emerges, which moves slowly in a straight line along its axis of symmetry. This dipolar structure takes on the same appearance as the dipoles observed in a linearly stratifed fl uid by van Heijst and Fl6r, 6 Voropayev et a/.,7 and

---• ----•--....., ••••••••••••• -• (a) (b) .... . (c) �- (d)

FvHF.

Typical examples of the vertical density profi le through a turbulent patch in a linearly stratifi ed fluid and in a two-layer fluid are presented in Figs. 4(a) and 4(b), respectively. In a linearly stratifi ed fl uid, the turbulent ed dies over the entire patch are competing with roughly the same restoring buoyancy force, which results in a single mixed region [see Fig. 4(a)]. For the two-layer case, the entrained fluid from the upper and lower layers have a buoyancy 8b/2 with respect to the turbulent eddies in the patch. The entrained lighter fl uid from the upper layer will be carried downward and the heavy fluid (entrained from the lower layer) will be carried upward by the overturning fluid of density ( p 1 + p 2 ) /2. Initially, the lighter fl uid tends to concentrate in the upper part of the patch, and the heavier fl uid will be in the low part of the patch [see Fig. 4(b)]. When eddies continue to overturn fluid across the midplane, they will experience an almost twofold increase of the buoyancy force, which may suppress the overturning if the inertia forces of the eddies are too weak to perform overturning motions against the increased buoyancy In some experiments it was observed that remnant stand ing interfacial wave motions interact with the dipole vor tices. In contrast to the traveling internal waves, these waves decay with a much longer time scale. The standing internal waves induce pulsations in vortex size, which cause stretching and compression of the vortices. This mechanism leads to changes in the vertical vorticity com ponent of the dipolar structure. As a consequence, the translation velocity and the horizontal size of the dipolar

0 .----.-------•-----•-------.
7.----------------. h(cm) l (a) 6 5 4 3 • 2 • 1 • • • • • • • • . . • • • • • • . . • • • • 0 �--�--�--�--�--�--� 0 2 3 4 5 6 t (s)
structure show oscillations around their mean values. Eventually, these standing waves are damped out and the dipole translates at a low velocity.

IV. COLLAPSE OF THE TURBULENT REGION

For a typical experiment, the maximum vertical patch size and the maximum horizontal spreading of the region were measured from video recordings. The vertical patch size as a function of time is shown in Fig. 5 (a), whereas the combination of the maximum horizontal and vertical patch sizes is presented in Fig. 5(b). Initially, the turbulent re gion grows equally in horizontal and vertical directions due to entrainment of ambient fluid. Once the infl uence of the buoyancy forces sets in, the patch growth is retarded and finally reaches its maximum size. Subsequently, the patch collapses and decreases in vertical size while the hor izontal dimensions increase continuously with approxi mately the same growth rate as before the collapse [see Fig. these, however, are considered to be of minor importance to the evolution of the vertical patch size, and are thus neglected in the scaling analysis, which will be discussed below.

25 .------------- em 20 15 10 5 • (b) Oo 0 0 0 0 • • •• • • •• ••• • • • • • •• 2 3 t (s)

A. The growth of the turbulent region

During the initial growth, the turbulent motions are dominated by the inertia forces. Thus, a scaling analysis similar to that advanced by FvHF for the growth of a turbulent region in a linearly stratified fluid can be used.

The growth of the turbulent region lasts for a very short period of time, so that the kinetic energy during that period can be assumed to be a constant, viz.,

a V a U6 u � v{i )' ( 2
)
where U0 is the injection velocity, V0 is the injected volume at time t=O, and u'2(t) and V(t) are the mean square horizontal turbulent velocity and volume of the turbulent region at time t, respectively. During the evolution, the turbulent region can be approximated by a self-similar con ical shape of characteristic dimension h(t), and thus its volume can be approximated by V(t) �h(t)3• For the tur bulence within the turbulent region it is assumed that u'2 w'2, with w'2 1 12 being the rms vertical velocity. Using the scaling for u'2 and V(t), one may write

12 V o U6 w =Bc�IT�, ( 3
)
where B 1

B. The maximum patch size

The growth of the turbulent blob continues until the vertical inertia forces (w'21h per unit mass) of the eddies have decreased to such an extent that the buoyancy forces ( f>b per unit mass) become of comparable magnitude. At that stage, the patch reaches its maximum vertical size, and the inertia forces are balanced by the buoyancy forces, viz., [START_REF] Van Heijst | Dipole formation and collisions in a stratifi ed fluid[END_REF] where he is the maximum vertical size of the patch. In other words, this relationship implies that the patch col lapse occurs when the bulk Richardson number he 8b/w'2 increases to a critical value B 4 . By using (3) and ( 6), it is possible to obtain

( V0 U6) 1 1 4 h c= B s � • ( 7
)
where B5 = ( B1 B4) 114• The experimental results for the maximum patch height are presented in Fig. 7, and it is obvious that they show good agreement with the proposed scaling. The best fit required for the data yields B5=0.5. By using (5) and (7) the time tc, at which the patch reaches its maximum size he, can be written as ----------------------------- I.

(

) 8
where B6=(B5/B3)512=1.7 is obtained by using B3 and B5• Direct measurements of the collapse time tc from the experiments yields a value of B 6 =2.2±0.6. This disparity can be attributed to the difficulty of locating the time tc to attain the maximum he. Nevertheless, the constants are of the same order, showing a reasonable agreement of the scaling arguments with the measurements.

C. The collapsing stage

After the patch has reached its maximum height, the eddies of the integral scale are too feeble to entrain the ambient fluid. An isolated patch of mixed fluid with a den sity different from the surrounding fluid remains, while the turbulence within the patch is weakened and fi nally is transformed into small-scale wave-like motions. Simple hy drostatic calculations show that there will be a horizontal pressure gradient between the inner and outer fluid, which causes the mixed patch to physically collapse and the mixed fluid to intrude in the ambient fluid at rest. The speed of the intrusion can be evaluated by a balance be tween the pressure gradient force and the inertia force of the intrusion: 16

where h(t) and r(t) are the total vertical patch size and the horizontal dimension, respectively, and C 1 is a constant. Assuming a double conical shape for the mixed patch, the horizontal spreading rate dr!dt can be estimated by h�=2?h(t). given in Table I.

By using ( 9), ( 10), and the expression for the maximum vertical patch size (7), it is possible to obtain 1 1

( 8b 7 ) 1 1 8 h-he = C3 v�ifo (t-tc), (11)
where C3=2 {iB5 3 12C 1• A comparison of data with (11) is presented in Fig. 8. The values of C3 for individual ex periments are listed in Table I. The constant c3 = 0.30 ± 0.09 is obtained by taking the mean value of c3 over all measurements. Although the data for each indi vidual experiment are in good agreement with the pro posed scaling, there is some spread of the values of C3 over all the experiments. This is probably due to the interaction of the turbulent patch with interfacial wave motions, and the neglect of the development of the density gradient within the patch during initial mixing. The value C 1 can be evaluated by using c3 and yields cl =0.04.

D. The Thorpe-scale measurements

In a turbulent stratified fl uid heavy fluid elements are displaced above the lighter ones by turbulent stirring mo tions, causing density inversions. One method to estimate the scale of such overturning motions is to determine the vertical displacements of fluid elements relative to a stable, monotonic density profile. This density profi le can be ob tained by rearranging the measured instantaneous profile into a stable one by displacing fluid particles to their neu tral buoyancy depth by a so-called "bubble sort" routine. 10 The required displacements are called Thorpe displace ments, which provide information about the overturning motions of the eddies. The rms value of these displace ments is known as the Thorpe length scale Lt. In order to measure the Thorpe scale of a decaying turbulent region, a number of experiments were performed under identical conditions, i.e., the same injection velocity and volume, and stratifi cation. For each experiment, two traversing conductivity probes-with a distance 2.9 em apart-were positioned at a given distance from the injec tion nozzle. These probes were shot vertically downward through the fluid at the instant that the turbulent region was at the same position as the probes. Information about the spatiotemporal evolution of the turbulent region was obtained by changing the position of the probes for each experiment. By using the measured density profi les, the maximum Thorpe displacement L1 max, the Thorpe scale L1, and the vertical patch size L P were determined: L P was defi ned as the region that has displacements exceeding 5% of L1 m ax. This threshold of 5% was necessary to suppress the noise peaks generated at the start and end position of the probe traverse. In addition, the patch size h was mea sured using video records, and the results are shown in Fig. ing the collapse of the patch the difference between these two patch size measurements becomes larger. Under these conditions, the color visualization technique indicates the physical size of the collapsing region, whereas the Thorpe displacement technique indicates the return of buoyant fl uid parcels to their respective stable position. The latter occurs faster, as was indicated by the measurements of FvHF.

The Thorpe scale L1 and the nondimensional Thorpe scale L/ Lrmax are plotted against time in Fig. 10. The graph shows that L1 decays approximately linearly with time. The normalized Thorpe scale L/ L1 max appears to • remain constant ( =0.4) during the evolution of the patch. This is in fair agreement with the previous observations of Itsweire [START_REF] Itsweire | Measurements of vertical overtes in a stably stratified turbulent fl ow[END_REF] and DeSilva and Fernando, 18 which showed Lrmaxl L1 = 2.75, for decaying stratifi ed turbulence and forced stratifi ed turbulence, respectively. Shown in Fig. 11 is a plot of L1 vs L P ; a linear increase can be seen with L1 Thus, the effective ratio is Ltmaxl Lp=0.45/ (Lp /2) =0.9.

This value is somewhat higher than that of DeSilva and

where z1 and oz are the level and thickness of the interface, respectively. When the wave amplitude is much smaller than the layer thickness, it is possible to obtain the disper sion relation as

k,= (2n+ l) a l!oz o}-8b!8z _ 4(n 2+ n) a l /8 � ]112 ol-Sb/8z ' ( 13
)
where n is the wave mode and k11 is the wave number (see, e.g., Ref. 21). Waves interacting with the dipolar vortex should have at least a wavelength of twice the diameter of the dipole in order to stretch and squeeze the structure as a whole; from Fig. 13, it follows that the dipoles have a mean diameter of approximately D=35 em, and thus waves with wave numbers k11<.2TT/2D=.9 m -I are the pos sible candidates. The frequency of these oscillations mea sured from Fig. 13 is approximately w=0.25 rad/s. Figure 14 shows that-according to relation ( 13 same frequency as the observed oscillation frequency (0.25 rad/s) of the dipolar vortex. Apparently, it is the wave mode 2 that causes the dipolar oscillation.

VI. CONCLUSIONS

A laboratory study on the evolution of a turbulent region injected at the interface of a two-layer fluid was performed. An isolated blob of turbulent fl uid in a two layer fl uid reaches its maximum height when the turbulent eddies in the blob come to a balance dictated by the vertical inertia forces of the eddies and the buoyancy forces asso ciated with them. Before this occurs, the effect of stratifi cation appears to be overshadowed by the inertia forces. After reaching this maximum height, the patch starts to collapse. This collapse leads to the formation of pancake like eddies, which may be modulated by interfacial waves. The results can be summarized as follows.

(i) The vertical thickness of the region grows indepen dent of the stratification for a time period tc= 1. 7 ( V0U 61ob 5 ) 11 8 • (ii) At tc, the stratification becomes important and arrests the vertical growth; the vertical size of the patch under this critical condition is determined by a balance between vertical inertia and buoyancy forces of the eddies, and is given by hc=.0.5 ( V0U 6/8b) 114• (iii) Measurements indicate an initial increase of the Thorpe scale to a maximum value (coinciding with he), and a decrease thereafter. The decrease is slow, but is still faster than the decrease of the patch size, owing to the physical collapse. The ratio of the Thorpe scale to the max imum Thorpe displacement was found to remain constant over the time period investigated. The relationship between the Thorpe scale and the patch size was measured to be L1::::;0.18L P " The maximum Thorpe displacement was found to be given by L1max1 L P = 0.45. This' is somewhat low compared to that observed in previous studies, and the development of stratification within the mixed layer during its evolution may be the reason for this observation.

(iv) The dipolar structure observed in the present two layer experiments has the same appearance as in a linearly stratified fluid. However, in a two-layer fluid internal wave motions appear to be more energetic than in a linearly stratified fluid because of their trapping at the interface; thus the interaction between the waves and the vortices is strong. These internal waves are mainly of mode 2, and subsequently stretch and squeeze the vortices, inducing os cillations of the translation speed and of the horizontal size of the dipole.

( v) The growth of a turbulent patch in two-layer and linearly stratifi ed fluids (studied by FvHF) appear to be governed by similar dynamics, but some differences can be observed in detail. For example, in linearly stratified fluids the entrained fluid into the patch is subjected to mixing over the entire turbulent fl uid column, whereas in the two layer case a sharp interface may develop at the plane of symmetry, thus inhibiting the exchange of fluid between upper and lower halves of the patch. Thus the overturning distances can be smaller in the latter case, and the con stants that appear in the different balances (scaling argu ments) can be different for the two cases.

  FIG. 1. A sketch of the experimental setup.

FIG. 2 .

 2 FIG. 2. A sequence of photographs, showing the evolution of the turbu lent region as viewed ejom the side. A volume of fluid V0 = 2 ml is injected horizontally at the level of the interface with U0=3.16 m/s and 8b=0.60 rn!s2• The pictures are taken at (a) 0.44, (b) 1.0, (c) 1.44, and (d) 44.8 s after the injection. The top-bottom distance in each picture represents the scale and is 30 em.

•FIG. 3 .

 3 FIG. 3 .. Soliton-like disturbance at the nose of the intrusion just before separation from the intrusion (a) at t=9.6 s, after separation and colli sion with the wall (b) at t=17.0 s, and at t=23.4 s. The experimental parameters for this experiment were 8b=0.615 m/s 2 , U0=3.6 m!s, and V0= 16.8 ml. (This particular experiment was carried out later in Eind hoven with a less powerful injection mechanism, which allowed only for smaller injection velocities. To obtain the same collapse height and there with a similar soliton, a larger volume had to be injected.) The dimen sions of the tank used for this visualization study were 90 X 40 X 40 em.The top-bottom distance in each picture represents the scale, which is 27 em.

FIG. 4 .

 4 FIG. 4. Density profile of the turbulent patch before the onset of collapse in {a) a continuous ambient stratification and in {b) a two-layer fluid. The stratifi cation parameters are the Brunt-VliisliHi frequency N=[(g!p)(8p/8z)]112=2.2 rad/s for (a) and 8b=0.48 m/s 2 for (b). The profiles were measured at (a) t=0.9 s and (b) 1.4 s after injection.

5 6 FIG. 5 .

 565 FIG. 5. Typical evolution of (a) the vertical patch size and (b) the corresponding growth in horizontal size compared to the vertical growth. The vertical patch size is denoted by the symbol 0 and the horizontal patch size by the symbol e. The experimental parameters fur this experiment were ob=0.757 m/s2, U0=5.44 m/s, and V0= 11.1 mi. The dimensions of the tank used for this visualization study were 300X40X40 em.

  5 (b)]. For all the experiments the vertical size of the patch is measured, until changes were no longer noticeable, i.e., when the horizontal size is approximately 15 em, and the influence of the tank width is still negligible. The collapse continues until the patch is confined in a thin layer of fluid, and pancake-like eddies start to form. In some cases it was observed that interfacial wave motions cause the mixed patch to oscillate with a period of about 0.5-1 s. These wave motions induce oscillations in the vertical patch size during the growth and subsequent collapse [see Fig.5(a)];

FIG. 6 .

 6 FIG. 6. The normalized vertical size of the turbulent region as measured for 16 experimental runs, with p= ( UOofV�) 1115• The experimental param eters and the proportionality constant B3 for each individual experiment are given in TableI.

FIG. 7 .

 7 FIG. 7. The maximum vertical size of the turbulent region as measured in 16 experimental runs. The experimental parameters are presented in TableI.

( 10 )FIG. 8 .

 108 FIG. 8. The experimental results for the collapse stage plotted in the form of Eq. (11), with q=(c5b71V6ifo>!!8• The experimental parameters are

9 . 6 FIG. 9 .

 969 FIG. 9. A comparison of the decay of the vertical patch size L P obtained from Thorpe displacements and the patch size h from video recordings, denoted by the symbols • and 0, respectively.

FIG. 10 .FIG. 11 .

 1011 FIG. 10. The Thorpe scale L, and the nondimensionalized Thorpe scale L,l L,max against time, denoted by the symbols • and 0, respectively. The data were obtained in an experiment with U0= 5.96 m/s, V0=4.5 ml, and eb=0.53 m!s1.

Fernando, 18

 18 Fig. 12 (a)]. By continuity of mass this diverging flow of mixed fluid induces a low pressure in the center of the collapsed region. After the collapse, it is observed that the intruded fluid flows back to this center of collapse. At this stage larger horizontal vortices are seen to arise, and their merger leads to the formation of the dipolar vortex struc ture. Simultaneously, the thickness of the dipole becomes larger than that of the flat collapsed region [see Figs. 12 (a) and 12(b)]. This dipolar vortex structure has an appear ance similar to the dipoles observed in linearly stratified fluids (see Refs.6 and 7).

6 FIG. 12 .FIG. 13 .

 61213 FIG. 12. Sideview photographs of the collapsed region (a) before and (b) after the emergence of the larger vortex structure. The photographs were taken at (a) t=9.0 s and (b) t=20.0 s, respectively. The experimental parameter values are 8b=0.615 m/s2, U0=3.6 m/s, and V0= 13.4 mi. The dimensions of the tank used for this visualization study were 90 X40X40 em. The top-bottom distance in each picture represents the scale, and is 28 em.

20 �FIG. 14 .

 2014 FIG.14. The plot of Eq. (13) for the modes n= 1, 2, 3, and 4.

TABLE I .

 I Table of the injection parameters and the evaluated constants of the laboratory experiments. U0, V0, and lib are the injection velocity, injection volume, and buoyancy jump across the interface, respectively. The experimental constants for the initial growth and the collapse B3 and C3, respectively (see the text) are given for each experiment, together with their correlation coefficients R2 of the approximations with the measurements.For experiments 20-22 a soliton-like wave motion occurred, and no regular collapse was observed.

	exp	U0 (m/s)	V0 (m1)	ẽb (m/s2)
	6	5.82	3	0.779
	8	5.82	5	0.779
	13	5.82	5	0.588
	15	2.92	7	0.653
	16	11.71	4	0.653
	18	9.08	5	0.653
	19 20 21 22 23	5.78 5.78 5.78 5.78 3.21	7 9 1 1 5	0.653 0.653 0.653 0.653 0.653
	24 25 26	10.81 1.8 6.4	5 3 3	0.653 0.371 0.371
	27	6.4	4	0.371

ACKNOWLEDGMENTS

We gratefully acknowledge the help of Dr. Prabath DeSilva for the use of software for density and Thorpe scale calculations, Mr. Leonard Montenegro for technical assistance with the conductivity probes and Professor J. C. R. Hunt for his critical remarks on various aspects of this paper. One of us (JBF) gratefully acknowledges fi nancial support by the Foundation for Fundamental Research on Matter (FOM) of the Netherlands Organization for Pure Research (NWO) and the U.S. Office of Naval Research (ONR), which made a seven-week visit to ASU possible.