CFD Analysis of the Particle and Melt Flow Behavior During Fabrication and Processing of TRIP-Matrix-Composites

Computational fluid dynamics simulations are conducted to supplement experimental investigations in order to gain a deeper understanding of physical effects during fluid flow and heat transfer in steel making and processing. This offers the possibility to examine physical effects of the liquid steel in greater detail and isolated of entire processes. Liquid steel is present in the fabrication processes of TRIP-Matrix-Composites, namely gas atomization of the steel melt to produce powder, and liquid steel infiltration into ceramic structures. It also occurs during further processing, e.g. in welding or coating. Numerical simulations of these processes are performed with the finite volume method using the free open-source software package OpenFOAM. The libraries are extended where needed. This includes formulations for phase change, heat sources, latent heat, additional forces, calculations for material properties in multiphase flows, and particle tracking. The models are used to simulate electron beam welding, infiltration, gas atomization, and flame spraying and to reveal significant effects for each particular process.


Introduction
Producing TRIP-steels and TRIP-Matrix-Composites through casting and powder metallurgy involves many physical effects both in the solid and molten state. This also applies for further processing like welding or flame spraying. A detailed understanding of these physical effects is crucial to achieve desired properties of the final product. Besides various experimental investigations, which are thoroughly described in this book, numerical simulations can be carried out to enhance this understanding.
The infiltration of ceramic structures like foams or bulk material is one route of producing TRIP-Matrix-Composites. In a second route, powder metallurgical processes, like e.g. atomization (Chap. 3), sintering (Chaps. 6 and 9) and powder forging (Chap. 8), are used. In order to produce the steel powder, gas atomization of the steel melt is performed. In both routes, the behavior of the molten steel is one main cause for the final material character. During further processing of the material, steel melts particularly occur during coating and joining processes. For the former, flame spraying can be used for metal matrix composites (MMC). For the latter, electron beam (EB) welding was found to show excellent weldability for TRIP-steels [1] and was therefore used.
For the numerical simulation of infiltration, atomization, flame spraying, and EB welding, computational fluid dynamics (CFD) simulations are used. The aim of these simulations is to create numerical models of the actual processes in order to reduce the amounts of costly experiments and to gain additional insights into the melt flow behavior and further related physical effects. In order to represent processes as effectively as possible, the right choice of fundamental equations to be solved and correct numerical model setups are crucial.
All CFD simulations of the described processes are based on the numerical solution of the conservation equations of mass, momentum and energy. The mass conservation equation for incompressible fluids with constant density has the form: ∇ · u = 0, (18.1) with the velocity vector u. The momentum conservation equation for incompressible fluids has the form ∂u ∂t + ∇ · (uu) = − 1 ∇ p + ∇ · (ν∇u), (18.2) where t is the time, is the density, p is the pressure and ν is the kinematic viscosity.
Since melt flow problems are accompanied by large temperatures and temperature gradients, the energy conservation equation is also solved to take the influence of temperature into account: with the temperature T , the heat capacity c p , and the thermal conductivity k. Equations (18.1)-(18.3) form a system of partial differential equations, which is solved by means of the finite volume method (FVM) [2,3] within the free, open source CFD software OpenFOAM [4].
If in addition to the steel melt, other fluid phases like gases are involved, the phases need to be distinguished. This is accomplished by using the volume of fluid (VOF) method and solving the transport equation for the phase fraction α: ∂α ∂t + u · ∇α = 0. (18.4) A very detailed description of the VOF method employed in OpenFOAM is given by Klostermann et al. [5].
For the particular processes, (18.1)- (18.3) were extended by source terms and other refinements. These extensions and their applications will be shown in the following sections for infiltration, gas atomization and EB welding.

Infiltration
Besides powder metallurgical processes, combination of steel and ceramic structures on the macro-scale is another possibility for the creation of TRIP-Matrix-Composites. In this production route, two different main approaches are investigated. One is the casting of liquid steel into molds that are filled with ceramic structures or particles. The other is the infiltration of ceramic foam structures with liquid steel (Chap. 4) [6]. The latter was investigated by means of CFD intensively by Klostermann et al. [7,8] and Schwarze et al. [9,10]. In order to perform beneficial CFD simulations, different steps need to be considered. These are efficient meshing strategies [8], steel flow and heat transfer inside the foam structure [7], and melt surface dynamics for the 2-phase-flow system [9,10].

Meshing Strategies
Sophisticated meshing strategies are the first step to efficient CFD simulations. The mesh is used to generate grid points in order to discretize the fundamental equations (18.1)- (18.3). The mesh offers a significant impact on the quality and timeconsumption of the simulation. In order to evaluate a meaningful meshing strategy for the discussed infiltration simulations, two different approaches were examined. A mesh, produced by a semi-automatic generation algorithm for hexahedral meshes, developed by Klostermann et al. [8], is compared with a mesh generated by a fullautomatic tetrahedral meshing technique. The comparison is performed for the flow in a Kelvin structure.
The main advantage of full-automatic techniques with tetrahedral or mixed tetrahedral/hexahedral grid cells is its independence of the geometry. Any forms can be meshed and automatically refined at critical areas. As an example, Fig. 18.1 shows a CT scan of realistic ceramic foam (Fig. 18.1a) and the respective mesh ( Fig. 18.1b). The automatic meshing also leads to quickly accessible meshing results. However, when the investigation of specific physical phenomena is desired, abstracted and well defined geometries are often examined. In these cases, structured meshes might lead to faster simulations and possibly more accurate results. To keep the meshing effort small, it is useful to employ automation for the geometry and mesh generation as much as possible. Since the Kelvin structure is a popular approach to approximate foam structures, Klostermann et al. [8] developed a strategy that does both geometry  [8] and b the automatically generated unstructured mesh with mesh refinement near foam boundaries [9] Fig. 18.2 A section of the generated Kelvin cell structure with triangular struts of one Kelvin cell (black) and one quarter of a Kelvin cell volume (dark grey) [8] and block-structured mesh generation for composed macroscopic forms that consist of single Kelvin structures. Figure 18.2 shows an excerpt of the geometry. For a detailed explanation of the meshing strategy, the reader is referred to Klostermann et al. [8]. semi-automatic block-structured mesh [8] In order to directly compare the two meshing approaches, the same geometry was meshed. Figure 18.3 shows sections of the two approaches. Both meshes were created with a similar number of grid cells. Steady, incompressible flow simulations through the Kelvin structure were performed to compare the mesh qualities, simulation errors and calculation time. Firstly, the block-structured hexahedral mesh mostly showed lower errors for different error estimates and two different Reynolds-number flows. Only in some cases, the errors of the unstructured tetrahedral mesh were slightly smaller. Secondly, the necessary number of iterations that was needed for convergence, was a lot smaller for the structured hexahedral mesh in all cases. Thirdly, the calculation time per iteration was lower with the hexahedral mesh. All in all, the performance of the hexahedral mesh was up to 94% better than of the tetrahedral mesh. In order to minimize computational costs the block-structured technique was consequently used for further studies of flow infiltration into Kelvin structures.

Mesoscale Flow in Kelvin Structure
Meshing the pores of the ceramic foam, as shown in Sect. 18.2.1 is very useful if only the melt flow is simulated. In order to represent the process of infiltration more physically, heat transfer between the liquid steel melt and the solid ceramic foam struts needs to be accounted for as well. To track the heating of the foam itself, the solid structure is meshed. Simulation domain with boundaries, regions and mesh for the fluid zone (light grey), primary meshed with hexahedric cells, and the solid structure (dark grey), primary meshed with tetrahedric cells [7] Therefore, the mesh was extended to the ZrO 2 struts in order to simulate complete heat transfer between the infiltrating liquid steel melt and the solid ceramic foam [7]. The aim of the performed CFD simulations was to investigate the flow and the heat transfer on the pore-scale (mesoscale). The single-phase flow focused on thermal effects with constant material properties and neglected the free surface flow and solidification. The equations for mass (18.1), momentum (18.2) and energy (18.3) were solved. The energy equation (18.3) is reassembled with the thermal diffusivity a = λ c p to: Additionally, a heat transfer equation for the solid region was solved: The following explicit boundary conditions were used at the interface between the melt and the struts: Temperature fields in the liquid zone (left) and the solid struts (right) after a t * = 0.49 and b t * = 6.14 melt through time [7] For details about the numerical model setup and the thermophysical properties, the reader is referred to Klostermann et al. [7].
The flow field is represented in Fig. 18.5. On the left hand side, the flow is visualized by streamlines, which illustrates the channeling of the flow in two alternating patterns, straight and meandering. Recirculation is visible through swirls behind the struts. On the right hand side, these recirculation areas are quantified, where the normalized velocity U * z in z-direction is negative, and therefore white. In realistic infiltration processes with a free melt surface, these recirculation areas could result in cavities as gas entrapment would possibly occur. Figure 18.6 shows the temperature field of (a) the flow field and (b) the struts for the melt through times t * = 0.49 ( Fig. 18.6a) and t * = 6.14 ( Fig. 18.6b). The ceramic was preheated to 1080 K. While the temperature field followed the fluid flow, the heating of the solid struts was delayed. Both behaviors were as expected. The temperature field of the liquid steel also revealed areas of lower temperatures in the recirculation areas behind the struts. We expect the melt to solidify in this area first, where gas entrapment is very likely as well.

Melt Surface Dynamics
In the previous considerations, the free melt surface and the gas phase have been neglected. An additional investigation of the multiphase system is, however, essential to cover further effects like cavity development. In order to investigate the behavior of the free surface of the steel melt in the porous ceramic foam during infiltration and to distinguish between the fluid phases, the model was extended by the VOF method [see (18.4)]. The material properties density and dynamic viscosity η were averaged between the steel melt phase (index m) and the gas phase (index g) depending on α in each grid cell [9]: (18.10) Figure 18.7 shows a simulation of the infiltration of liquid steel into the ceramic matrix with a highly resolved surface tracking [9]. The gas is noticeable as dark entrapped bubbles in between the lighter steel melt. The bubbles form as a consequence of the rapid infiltration behind the struts.
The detailed simulations came with tremendous computational effort because very small time steps t were necessary. The required time step for free surface flow is calculated by the following condition [11]: (18.11) where C η is a constant and x is the cell length. As described by Schwarze and Klostermann [10], this time step condition is much more dominant than the usual Courant-Friedrichs-Lewy-condition, which reads: Here Co is the Courant number, which is in most cases Co ≤ 1 to ensure stability, and u is the velocity magnitude. In (18.11) C η most widely depends on numerical treatment of the model equations and therefore the time step is very sensitive to adequate implementations. Schwarze and Klostermann [10] have found that the material parameter averaging between the phases have a particular impact. Thus, an improved strategy was found for the viscosity interpolation and the harmonic mean was adopted: Testing different time steps with the different viscosity averaging equations reveals their impact on the stability of the simulation. In Fig. 18.8, the standard (18.10) and the improved (18.13) viscosity averaging algorithms are compared. For the low timestep t 1 both implementations worked stable. By increasing the time step to t 2 , the standard implementation destabilized. However, with the improved viscosity   Comparison of the free surface between simulation (left) and model experiment (right) for two different through times a t * = 0.199 and b t * = 0.364 [10] averaging, the results were comparable to the smaller time step and the simulations remained stable. Figure 18.9 shows a qualitative comparison of the infiltration process between simulation and experiment, It reveals a very good agreement.

Atomization
Alongside infiltration of macroscopic structures, powder metallurgy is the main route for the fabrication of TRIP-Matrix-Composites. The controlled production of steel powder with desired properties and particle size distributions is one of the major challenges. In the light of this, CFD models were developed and simulations were performed to gain an understanding of parameters that influence the development of steel droplets during the primary breakup of molten steel [12]. In order to apply additional models like secondary breakup or applications like flame spraying, particle size distributions need to be tracked. By converting the continuous Euler phase of the VOF simulations to disperse Lagrangian particles, this can be accomplished. For that purpose, a particle tracking and conversion model has been developed and implemented [13]. The application of tracked particle size distributions for powder flame spraying is shown accordingly.

Influence of Process Parameters on Primary Breackup
The production of MMC using a sintering process needs analog particle sizes for the ceramic and steel powders to obtain a high green strength [14]. Atomization, the predominant production process for steel powder, subdivides into two main tech-niques, the free fall design and the close coupled (confined) design [15]. The main difference of the designs is the positioning of the nozzles, which deliver the inert gas for the atomization of the steel melt jet. Therefore, different atomization states are feasible. In case of the closed coupled design, a prefilming process of the metal melt jet develops due to the flow structure of the inert gas, which finally atomizes the melt into droplets. On contrary, a compact metal melt jet develops for the free fall design, which is atomized by an inert gas stream at a defined distance. Due to the inert gas atomization, steel powder particles with diameters in range of 0.04-2000 µm are generated [16]. A method to enforce the atomization process is the direct manipulation of the free falling melt jet. Therefore, a primary breakup process is forced by defined disturbances for the inflow of the melt in the nozzle. This leads to a pre-disintegration of the metal melt jet into ligaments and larger droplets [12]. In general, an experimental investigation of the steel melt atomization is very complex. Mainly, the high melt temperature limits the measurements [17]. The management and control of the process of metal melt gas atomization is a kind of 'black art' due to the demanding process conditions and requirements [18]. Furthermore, experimental setups require a certain durability and stability as well as optical accessibility for a detailed analysis of the process. Consequently, numerical simulations are a meaningful alternative to research metal melt atomization processes to shed light on specific influencing attributes.
Critical factors for primary breakup of the metal melt jet using a free fall design are constructive parameters (e.g. the nozzle design, nozzle obstruction via the stopper rod) and operation conditions (e.g. rotation of the melt flow, nozzle inlet flow disturbances). Basically, a numerical investigation determines their influence to the metal melt jet behavior downstream of the nozzle. Figure 18.10 shows sketches of a tube nozzle ( Fig. 18.10a) and a divergent-convergent industrial nozzle ( Fig. 18.10b). The presented numerical simulation depends on different assumptions such as: • three-dimensional and transient flow behavior, • fully developed laminar melt inflow, • the flow regime is in transition range between laminar and turbulent flow, • isothermal and incompressible conditions with constant material properties. Therefore, the governing equations for continuity (18.1) and momentum (18.2) are solved. The turbulence is not modeled with an explicit subgrid model. Instead, implicit large eddy simulations (ILES) were performed, which generate the dissipation of the small structures by a specific spatial numerical discretization scheme of the convection terms in the momentum equation. By the means of the VOF, introduced in Sect. 18.2.3 [see (18.4), (18.9), (18.10)] the dynamics of the phase interface between the melt and gas is tracked during the simulation. Figure 18.11 shows the a simulation result of a temporal metal melt jet behavior downstream the industrial nozzle. The nozzle outlet diameter is consistent D = 4 mm for both geometries (see Fig. 18.10). The axial nozzle inlet velocity is constant U ax = 1.55 ms −1 based on the fluid column of the melt inside the tundish. Additionally, a circumferential velocity of u rot = 1.83 ms −1 is defined, which mimics the rotation of the melt inside the tundish caused by the magnetic heating process. For detailed information of the specifications, the reader is kindly referred to the work of Neumann et al. [12]. In Fig. 18.11, the contour plot of the jet is depicted by the phase fraction α = 0.5 for the industrial nozzle. It should be noted that the interface between the melt and gas phase smears over 4-5 cells [8]. However, respecting this criterion a high computational grid resolution is defined. Additionally, the gas velocity field is shown in the background. In general, the temporal behavior of the melt jet shows the primary breakup of ligaments and larger droplets. A helical structure superposed by a hollow-cone jet shape develops, whereas close to the nozzle a compact almost cylindrical jet occurs with less surface waves. Here, the disintegration is directly driven by the circumferential velocity of the steel melt, which shows a high potential to improve the secondary atomization process due to the high excitation of the melt jet. Therefore, Fig. 18.12 displays a detailed analysis of the different nozzle designs considering different circumferential velocities. Depending on the nozzle design different stages of the jet behavior can be observed. The tube nozzle ( Fig. 18.12a) shows for: • u rot = 0 less disturbances of the melt jet, which lead to surface waves and diameter changes, • u rot < u ax a stabilization of the jet leading to a smoothing of the secondary structures, • u rot ≈ u ax a helix shaped metal melt jet with less excitations onto the jet surface, • u rot > u ax an irregular jet behavior, which splits into two parts.
On contrary, the results for the industrial nozzle design (Fig. 18.12b) show a higher excitation level of the melt surface. However, there are also different jet behavior stages considering the circumferential velocity contribution. In case of u rot < 0.5u ax , the general jet behavior is preserved compared to the undisturbed variant. Surface waves, diameter changes and partial disintegration occurs. A helical jet shape develops by increasing the circumferential velocity to u rot > 0.5u ax . Furthermore, the jet tip moves with time and ligaments and larger droplets disintegrate. However, the surfaces waves are persistent. For u rot > u ax , the helical jet is superposed by a hollow-cone shape and the jet oscillates frequently. Higher velocity increases the spreading and separation of the jet. The ligaments and particles are in magnitude of the secondary structures. In summary, the nozzle design mainly influences the jet excitation due to surface waves. However, the decay of the melt jet is primarily forced by specific operation conditions. Deeper insights in this research can be found in the work of Neumann et al. [12].

Particle Tracking and Conversion
Originally, the VOF method is capable of resolving arbitrarily formed interfaces (e.g. water surface). However, a direct tracking of specific parameter (e.g. particle size) is not immediately available for VOF simulations. This is comparable to experimental research, where a specific evaluation process is required to identify the desired variable [19]. The same applies for VOF simulations, where a specific tracking routine manages this task. In case of steel melt atomization, droplet or particle related characteristics are of particular interest (e.g. particle size distribution). However, a VOF simulation of a steel melt atomization needs a sufficient mesh resolution to resolve the droplets, which are in the magnitude of 0.1-1000 µm [16]. This in combination with a multiplicity of droplets in the dense spray region would cause excessive computational costs [20,21]. Here, VOF simulations are limited to the complete primary breakup process, where the droplet sizes are reasonable for VOF. Secondary breakup processes are mostly described by an Eulerian-Lagrangian approach, where the com- Fig. 18.13 a Spatial CCL: Particles (grey) represented by the phase fraction field (left), particle face identification via phase fraction threshold (middle), grouping of cells by CCL to particle slices (right) and b temporal CCL: Particle movement through layer evaluation over time putational costs are much lower [22,23]. Therefore, the VOF droplet information (e.g. particle size, position and velocity) have to be transferred from the Eulerian to the Lagrangian simulation. Here, a capable tracking algorithm is indispensable. In general, there are two fundamental methods, namely three-or two-dimensional tracking [24][25][26][27][28][29][30]. 3D-tracking detects particles in the entire domain, while 2Dtracking defines a fixed positioned evaluation layer for the detection of the Eulerian particles. Both methods employ the connected component labeling (CCL) approach, which is based on the graph theory. However, whereas the 3D-tracking demands exclusively spatial CCL, 2D-tracking needs spatial and additionally temporal CCL for the Eulerian particle detection. In Fig. 18.13, the variants are outlined exemplarily. The basis for the tracking algorithm is the phase fraction α field given by the VOF simulation. Regarding the spatial detection ( Fig. 18.13a), particles have to be identified first. Therefore, an interrogation sequence of the phase fraction threshold has to be satisfied in order to be considered for the actual CCL step for the condition α > α t . Afterwards, neighborhood relationships between the detected cells have to be analyzed. Finally, individual particles are classified. The 3D spatial CCL directly detects an entire particle, whereas the 2D spatial CCL identifies slices of the particles due to the 2D evaluation layer (arbitrarily placeable). Therefore, an additional temporal CCL is important to connect multiple slices to one particle.
The simulations considering VOF metal melt atomization are performed in Open-FOAM v1812 [31]. Therefore, the officially implemented tracking algorithm, called the extract Eulerian Particles function object, is advanced to incorporate sophisticated source code for the temporal tracking. Figure 18.13b shows the correct temporal connection of the individual slices over ensuring correct neighborhood relationships within different particles.
To study the performance of the enhanced tracking algorithm, an extensive numerical validation of a well-known test case is intended where almost all circumstances   [13]. b Sensitivity analysis of the phase fraction threshold α t (Reference α t = 0.05, C d = 16 and Co = 0.1) of the simulation are controllable. Therefore, the test case of a moving water droplet in air is chosen using different setups, see Fig. 18.14. Different particle shapes (spherical, H-shaped) and number of particles (single, multiple particles) are analyzed defining different diameters d 0 and velocities u 0 . Furthermore, the influence of the mass and momentum equations are analyzed. Using the frozenFlow option of OpenFOAM enables simulations with constant domain velocity and negligence of the surface tension of the droplet. The simulation domain is defined by cyclic boundary conditions at the inlet and outlet, such that the particle can directly reenter the domain after leaving. This allows an iterative analysis of the tracking algorithm for constant particle conditions. In total, 30 cycles for each setup are performed. The mesh resolution is 32 cells per diameter. Figure 18.15a shows the simulation results of the VOF simulation with the enhanced 2D-tracking algorithm, with the normalized droplet volume V * and normalized displacement x * , which are described by: 14) within the box plot diagram. The diagram reveals a very good accuracy of the considered particle attributes. The discrepancy of all regarded setups is less than 3% for the normalized volume. Noticeable aberrations are existent for H-shaped droplet B) and the setup with multiple droplets of different shapes and sizes D). The normalized position of the particle is generally almost zero. However, if surface tension is considered (marked with *), the shape of the particle cannot be preserved and partially the droplet breaks up and shift in position. This leads to the higher discrepancies in case C) of Fig. 18.15a. Altogether, the enhanced tracking algorithm shows satisfying results for all tested parameters.
Since the reliability of the enhanced tacking algorithm is proven, a sensitivity analysis of the most influencing parameter should declare its limitations. Therefore, the influence of the phase fraction threshold α t , the droplet resolution C d = d 0 x and the Courant number [Co, see (18.12)] to the normalized volume V * are investigated and illustrated in Fig. 18.15b. The reference configuration (blue) is defined by α t = 0.05, C d = 16 and Co = 0.1. Alternative 1 (green) has a reduced droplet resolution of C d = 8 and alternative 2 (red) has a increased Courant number of Co = 0.5 compared to the reference. The remaining variables are equal to the reference. Figure 18.15b shows comparable trends for the different configurations. With rising phase fraction threshold increases the deviation of the tracked particle volume to the initialized. Regarding the reference configuration, the normalized volume differs ≈8% for the maximal tested phase fraction threshold. However, compared to the alternatives it shows less aberration and a low standard deviation. Less stringent adjustments for the droplet resolution and Co (green and red line) highly impair the particle volume detection. With rising phase fraction threshold, the difference to the reference results increases. Furthermore, the standard deviation is significantly influenced. In summary, the phase fraction threshold is the most critical parameter, which can highly influence the tracking results in VOF simulations. For further details, the reader is referred to the work of Spitzenberger et al. [13].

Flame Spraying
Flame Spraying is a coating technique, where steel and ceramic powders or sintered rods are molten to droplets and sprayed on a substrate. During flame spraying, the hot gas mixture is combusted in order to melt the coating material, which is either feeded Fig. 18.16 a Sketch of the nozzle with an annular air inlet and orifices for the oxygen-fuel mixture inlet and b front view of the 3D nozzle geometry as powder or a sintered rod. The air stream then atomizes the molten material and carries the droplets towards the substrate. For the basic principle of thermal spraying techniques, the reader is referred to the literature [32]. This technique is shown as an example for the need of Lagrangian particle tracking.
In the CFD model, the process is subdivided into two decoupled calculations for gas and droplet injection. At first, an oxygen-fuel mixture and air are injected through different inlets in the nozzle. An annular orifice is used for the air while the oxygenfuel mixture enters through twelve circular orifices. Figure 18.16a shows a sketch of the nozzle geometry, and Fig. 18.16b a 3D view of the nozzle used in the simulation. In the CFD model, combustion and melting are not taken into account. The hot oxygen-fuel mixture and cold air are injected in the nozzle and subsequently form a free jet. The hot mixture streams towards the substrate, whereas the temperature and velocity fields are calculated in a steady state simulation. An example of the resulting velocity field is shown in Fig. 18.17.
Since combustion, melting and atomization are not calculated in this model, Lagrangian particles are injected into the domain at the rod tip in the second part of the simulation. Therefore, the particle size distribution of the introduced melt droplets needs to be predefined. It can be estimated through upstream calculations of the gas atomization (Sect. 18.3.1) and tracking of the particle sizes (Sect. 18.3.2). This approach is most accurate for powder flame spraying processes, where the powder is used directly instead of the intermediate sintering into rigid rods. Therefore, in the case of rod flame spraying, Gaussian distribution with an estimated mean diameter can be used alternatively.
After introducing the droplets into the domain, they are carried by the precalculated gas jet toward the substrate and impact and stick to its surface. Figure 18.18 visualizes this process. The impact of the Lagrangian particles is only sufficient to estimate the location of the coating. In order to calculate the actual coating buildup, the Lagrangian point particles need to be transformed into particles with finite sizes that are able to deform and solidify after the surface impact. This procedure is part of ongoing studies. Furthermore, more accurate particle size distributions could be achieved by directly measuring the particles after flame spraying.  [33] and dissimilar joints [34,35]. The numerical models described here are based on their experimental work.
The basis for numerical investigations of EB welding are, on the one hand, underlying physical effects, and on the other hand, the experimental conditions. The former are described by fundamental mathematical equations and their translation into numerical equations via discretization on finite grids that reproduce the experimental domain. The latter are described by the structure and components of the domain and boundary and initial conditions. Recently, a comprehensive CFD study of the melt flow and keyhole dynamics was carried out by Huang et al. [36] with constant thermophysical properties. Furthermore, Borrmann et al. investigated the temperature sensitivity of the thermophysical properties of X3CrMnNi 16-7-6 TRIP-steel during EB welding [37] and the coupling between steel melt and ceramic particles in EB welding of TRIP-Matrix-Composites [38] in their CFD simulations.
The basic domain for EB welding simulations can be seen in Fig. 18.19. It sketches the essential parts of the welding process that need to be taken into account.
The welding seam is produced by heating the base material up to the melting point. As a result, a tear-drop-shaped weld pool forms. When the melt cools down at the trailing edge of the weld pool, it solidifies and forms the fusion zone. The in-depth weld pool and fusion zone shapes are illustrated in Fig. 18.20. It is evident that the geometry is highly anisotropic and therefore its three-dimensional character needs to be considered.
In the experiments, the weld pool is formed in consequence of the EB energy input, which has a highly transient behavior. In order to reduce computational costs to an acceptable level, the EB is not resolved as such. Instead, a heat source is  . 18.21 Sketch of the heat source (light blue) inside the weld pool (red). It is used for modeling the heat input, rather than calculating the electron beam itself operated, which mimics the EB heat input into the material. Figure 18.21 sketches its basic idea. The value and form of the heat source is constant throughout the whole simulation. However, the form of the weld pool is varying over time.
Detailed explanations of the physical effects and their numerical implementations are described in the following sections.

Phase Change and Heat Source Model
The increased complexity of the physical problem that arises as a consequence of heating and phase change requires additional source terms. For this reason, the fundamental equations (18.2) and (18.3) were extended, based on the elaborations of Voller and Prakash [39], Rösler and Brüggemann [40], and Miehe [41]. The momentum equation (18.2) is extended by the source term S m , which accounts for the solid fraction and damps flow in solid and partially solid regions: (18.16) and the source term S b , which accounts for the Boussinesq approximation, with which the density dependents linearly on the temperature in the gravity term: . (18.17) In the preceding source terms, d is a characteristic length scale of the solid material like e.g. a grain size, f s represents the solid fraction, which lies in between 0 and 1, c 2 is a numerical constant to prevent division by zero for f s = 1, g is the gravitational acceleration, β is the thermal expansion coefficient and T and T liq are the temperature of the fluid element and the liquidus temperature, respectively. Extending (18.2) with (18.16) and (18.17) leads to the final momentum equation that is solved in the EB welding model: To account for the latent heat of fusion H during melting and solidification, (18.3) is extended by the source term S , which is denoted as: (18.19) Additionally, the heat source term S q , which is sketched in Fig. 18.21 is added to the energy conservation equation: Combining (18.3), (18.19) and (18.20) leads to the form: (18.21) which is solved in the numerical model for EB welding. The heat source model is responsible for forming the weld pool during melting and subsequently the fusion zone after solidification. Because of the complex shape of the welding seam, which the heat source needs to mimic, it consists of a double ellipsoidal part for the weld nail head and a conical part for the middle and lower region, as proposed by Lundbäck and Runnemalm [42]. The double ellipsoidal part q e , recommended by Goldak et al. [43], was used: The conical part q c is: (18.23) Here, Q in is the heat, which is produced by the electron beam, a c/e , b c/e , and c c/e are the width, height and length of the heat source parts conical c and double ellipsoidal e, respectively. x, y and z are the Cartesian coordinates and d c is the narrowing of  the conical part, if its diameter is decreasing with depth. Figure 18.22 illustrates the heat source geometry and the related parameters. Combining (18.22) and (18.23) leads to with the weighting factor γ . Then, (18.24) is inserted into (18.20). The fundamental equations (18.18) and (18.21) were discretized with FVM with second order accuracy schemes. For thorough details about the numerical model, and the values of all thermophysical properties, the reader is referred to Borrmann et al. [37].
The introduced model was used to calculate the formation of the weld pool and subsequently the fusion zone of a welding seam that has been produced experimentally beforehand. The chemical composition of the welded base metal is shown in Table 18.1. Based on this composition, temperature dependent thermophysical properties were derived and transferred to the numerical model. Details about these properties and their importance are described by Borrmann et al. [37].
In order to determine target values for the outcome of the numerical model, welding seams have been produced within a vacuum chamber. The related parameters, which are deduced from the evaluation of the welds, are listed in Table 18.2. A cross section of the welding seam that is used here is shown in Fig. 18.23. The typical fusion zone geometry with upper nail head and long and narrow middle and bottom part is clearly recognizable. The weld reinforcement and root are also apparent in the experimental cross section. However, in the numerical simulation those are not considered.
Due to melt convection and heat conduction, the weld pool enlarges over time and therefore, the related heat source needs to be significantly smaller than the weld pool. Appropriate parameters are shown in Table 18.3. Figure 18.24 illustrates the coherence between heat source and weld pool. Due to the choice of heat source parameters, the upper double-ellipsoidal part of the heat source appears as a hemisphere, which is significantly smaller than the nail head of the weld pool. This applies for length, width and depth, likewise. The conical part of the heat source is noticeably smaller as well.
During the transient simulation, the weld pool and fusion zone are created gradually. Figure 18.25 shows process snapshots for four consecutive points in time. At the beginning of the simulation, the base material is heated exactly at the location of the heat source. The melting follows this behavior and therefore, the weld pool is shaped according to the heat source geometry in Fig. 18.25a. In the process of time, the weld pool widens and prolongs due to heat transfer through convection in the weld pool, heat conduction and further melting of the base material. This can be seen in Fig. 18.25b, c. Furthermore, the weld pool clearly solidifies and builds the fusion zone, which is indicated as a new area in blue. At later time steps, the geometry of the weld pool and fusion zone sustains, as seen in Fig. 18.25d. The cross section of the welding seam can be obtained once the weld pool is solidified and the fusion zone is formed. Comparing it to the experimental weld geometry, see Fig. 18.26, reveals the accordance between experiment and numerical model. Further examples for weld geometries can be found in Borrmann et al. [37].

Influence of Keyhole on Fluid Flow
In EB welding, the formation of a vapor column, the so called keyhole, plays an important roll for the formation of deep welds. The necessity to perform EB welding under very low ambient pressure, like in our experiments p A = 2 Pa, leads to a physical particularity during vaporization. The saturated vapor pressure p s,0 , which is the balanced pressure that determines the related vaporization temperature T s,0 , is  [37] significantly lower than the pressure p 0 with vaporization temperature T 0 under standard conditions. p s,0 is composed of the ambient pressure, the hydrostatic pressure of the steel column, and the dynamic pressure, which is affected by the melt flow velocity. Depending on the locally balanced p s,0 , T s,0 can be calculated by using the simplified Clausius-Clapeyron equation [44]: Solving (18.25) for T s,0 leads to: Since T s,0 is considerably lower than T 0 , vaporization occurs a lot earlier and closer to the melting temperature.
Due to the low ambient pressure, the hydrostatic pressure of the steel melt column gains in importance. Therefore, the vaporization temperature increases noticeably with increasing depth, due to a rising steel melt column. Additionally, the melting occurs very rapidly during EB welding with very high beam intensities. This leads to an overheating of the steel regarding its vaporization temperature T s,0 to an increased vaporization temperature T s . As a consequence of the overheating, the vapor pressure after vaporization increases noticeably. The now occurring pressure difference results in the so-called recoil pressure p r . It might as well be solved with the Clausius-Clapeyron equation: In dependence of the overheating temperature, the recoil pressure has a significant impact on the surface velocity of the melt around the keyhole. Figure 18.27 clarifies this dependence by showing the maximum keyhole surface velocity, depending on the overheating temperature T . A nearly linear correlation arises. Two example velocity profiles are given for T = 10 K (red circle) and T = 250 K (blue circle). The keyhole is indicated as the grey area within the weld pool. While the maximum weld velocity occurs at the keyhole surface in both cases, it is more than 20 times greater for the higher overheating temperature. Two additional effects are noticeable. Firstly, with increasing surface overheating, the keyhole tends to be tighter, because vaporization occurs at higher temperatures, whereby the keyhole has less time to develop. Secondly, the melt flow shows increased vorticity due to higher velocities and, as a consequence thereof, a higher Reynolds number. Since the simulations were performed without turbulence modeling, only vortexes that are resolvable through the mesh resolution are calculated.
The consideration of vaporization and the keyhole emergence, lead to an altered geometry of the weld cross-section in the simulation due to a widening of the weld pool. To correct this deviation, the heat source model parameters can be adjusted. In order to further refine the model, the front and rear length of the heat source,  . 18.28 Comparison of the solidified weld pool length in a the experiment with the top view of the simulated weld pool, b without keyhole and c with keyhole [37] indicated by f and r , respectively, were adjusted according to recommendations by Goldak et al. [43]. While the front part stays half the size of the weld width, the rear part parameter is increased to the double weld width. Nevertheless, the melt pool doesn't follow the heat source specifications precisely, because of further heat conduction, convection in the melt pool, and the recoil pressure. The refined heat source parameters are shown in Table 18.4. In order to restrict the infinite expansion of the Gaussian distribution, values with less than 5% of the integral size were truncated. The clipped amount was added to the integral distribution size accordingly. The weld pool length is indicated in Fig. 18.28. Figure 18.28a shows the experimental size, after solidification. The visible length of the solidified weld pool is 6.5 mm. Figure 18.28b, c present the weld pool sizes without and with the keyhole, respectively. While the welding seam cross section matches in both cases, the weld pool length agrees better to the experiment with using the keyhole and further heat source refinements. The relative error in case (b) is δ = 32.8% and in case (c) δ = 17.2%. The remaining error indicates that further improvements like free surface representation together with surface tension and corresponding effects need to be considered for more precise studies.

Dissimilar Welding of MMC-Steel
The weldabilty of TRIP-steel has been shown to work very well [1]. The situation is different, however, for the TRIP-Matrix-Composite. Steel and ZrO 2 powder is combined to form the MMC. The presence of ZrO 2 particles in the melt leads to effects that impair the weld result, like melt ejection and cavities. We assume that direct interaction between the particles and the electron beam and the interaction with the overheated melt pool and keyhole can lead to large energy input and potentially Fig. 18.29 Sketch of the dissimilar welding strategy with beam offset [38] even rapid vaporization of the ceramic particles. In order to find solutions to ensure the weldabilty of the MMC, experimental studies were performed [34,35]. Very good weldability was achieved, when the MMC plate was welded together with a particle-free TRIP-steel plate. By using a beam offset towards the latter, the electron beam as well as the melt pool are shifted. In that way, the disturbance of the ceramic particles is avoided. The main principle of this beam offset welding strategy is shown in Fig. 18.29. The weld pool cross-section is illustrated by the typical geometry in red. By shifting the beam as far as indicated, neither the beam, nor the emerging keyhole are interacting with the ceramic particles, illustrated as black spheres. However, they might still be present in the weld pool. In that way, they are separated from the steel matrix in the MMC and might float up or through the weld pool toward the particle-free steel. In order to investigate this behavior, the numerical model was further extended with the coupling between melt flow and ceramic particles. The necessary extensions are presented, together with parameter studies and results. For the complete numerical model and investigations, the interested reader is referred to Borrmann et al. [38].
The particle movement was driven by forces that are exerted to the particles by the melt flow. This was considered by solving their equation of motion, which reads Here, m p is the particle mass and u p is the particle velocity. The gravity force F G and buoyancy force F B are calculated with with the particle density p . The drag force F D is calculated by

Fig. 18.30
Top view of the calculation domain a before and b after welding simulation with base material (grey) and welding seam (red). One welding partner is filled with particles (blue), the other is particle-free [37] where C D is the drag coefficient, and A is the representative particle area. The particles were moved in each time step subsequently to the melt flow calculation. C D is defined with: Re ≤ 1000 0.424 Re Re > 1000. (18.31) Re is the particle Reynolds number, which is: After the field calculations for the melt are finished, particle forces are calculated. The particles are displaced at the beginning of each time step. Figure 18.30 shows the particle displacement. Figure 18.30a shows the top view of an initial particle distribution. The particles are generated at one half of the continuous domain on an equidistant mesh before the simulation starts. During the formation of the welding seam, particles were moved and displaced. Figure 18.30b shows the particle distribution after the end of the simulation. The welding seam, which connects the particle-free steel and the MMC, is now present. An offset of the seam reveals the offset of the heat source toward the particle-free steel. It is evident that ceramic particles were moved to a region inside the seam, where no particles were present previously.
Depending on the offset-ratio, it can be estimated if the particles would move toward the keyhole and interact with the steel vapor and the electron beam, which would still cause bad weldabilty. Therefore, different offset ratios were investigated. In Fig. 18.31, the total number of particles n p,d that were displaced and enriched in the previously particle-free steel is shown. Depending on the beam offset, the particle count differs significantly. The beam offset values are referring to the half of the nail head width such that for 100% offset, the welding seam is located completely on the particle-free welding partner. Between no beam offset (0%) and a small beam offset Fig. 18.31 Total number of particles that are enriched in the particle-free steel after welding in dependence of the relative beam offset [38] (25%), the total number of enriched particles doubles. Even though less particles are extracted from the MMC in the second case, the separation is higher than in the first case, because the particles were transported towards a greater region in the particle-free steel. For higher offsets, the trend reverses and the particle enrichment decreases as the count of separated particles is much smaller than for low beam offsets. One reason for this is the fact that for high beam offsets, the total weld pool area decreases disproportionately, because the offset has a higher impact on the displacement of the narrow area of the weld pool. For the highest beam offset (75%), the narrow welding seam part is completely located at the particle-free steel welding partner, which leads to a kind of braze welding. Therefore, the total number of displaced particles decreases largely.
The number of particles in the core area of the weld pool is of special interest. The core area is defined as a cylinder with half the diameter of the narrow part of the welding seam. A comparison between the state before and after welding is shown in Fig. 18.32. It is evident that the particle-count in the core area nearly remains the same for 0 and 10% beam offset. Increasing the offset leads to a shift toward a region, where no particles were present in the initial state. However, for 25 and 50% beam offset, particles were enriched in this area, after welding. While this is possible in the simulation, for real EB welding applications, particle would interact with the electron beam and disturb the process. Because of the particle enrichment, this would even be the case for the higher beam offsets. Only for the highest investigated beam offset of 75% no particles occurred in the core of the weld pool. This means that optimal welding conditions occur for high beam offsets, where the majority of the weld pool lies on the particle-free steel side. Because of the braze welding effect, the bond of the welding partners is still given [34].

Conclusion
The relevant processes for the production and further processing of TRIP-steel and TRIP-Matrix-Composites were investigated using CFD simulations. Various numerical models were implemented and existing solvers were extended.
Implementations into OpenFOAM have been made to improve stability, reduce computational costs in VOF calculations, to transform between different approaches, namely Eulerian and Lagrangian, and to calculate and couple particular physical effects. All implementations have been applied for chosen examples. The ability to calculate the infiltration of steel melt into macroscopic foam structures and free surface tracking to predict the possibility of infiltration and cavities has been achieved. The atomization model is able to predict primary breakup and reinforcing process parameters and track droplet sizes, position and velocity after primary breakup. Lagrangian particle distributions can be used to investigate sizes of coatings in flame spraying and the impacting positions and speed of droplets. The EB welding model is able to calculate weld pool shapes and melt velocities. Furthermore, the movement, depletion and accumulation of ceramic particles in the welded components can be evaluated. Through the calculation of the keyhole and therefore the penetration depth, the welding model is also able to estimate the weldability for varying EB parameters like beam intensity and feed rate.
Further extensions of the numerical models are possible. This includes the secondary breakup of the Lagrangian phase particles during gas atomization to get final particle size distributions. Furthermore, the coating buildup during flame spraying can be achieved with a reverse particle tracking algorithm which creates VOF droplets from Lagrangian particles. Models for solidification and surface tracking can be combined to include further effects like surface tension and Marangoni convection into the welding model and to evaluate their influence on particle-melt-interaction.
Open Access This chapter is licensed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license and indicate if changes were made.
The images or other third party material in this chapter are included in the chapter's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the chapter's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.