A simulation-based approach to characterise melt-pool oscillations during gas tungsten arc welding

Development, optimisation and qualification of welding and additive manufacturing procedures to date have largely been undertaken on an experimental trial and error basis, which imposes significant costs. Numerical simulations are acknowledged as a promising alternative to experiments, and can improve the understanding of the complex process behaviour. In the present work, we propose a simulation-based approach to study and characterise molten metal melt pool oscillatory behaviour during arc welding. We implement a high-fidelity three-dimensional model based on the finite-volume method that takes into account the effects of surface deformation on arc power-density and force distributions. These factors are often neglected in numerical simulations of welding and additive manufacturing. Utilising this model, we predict complex molten metal flow in melt pools and associated melt-pool surface oscillations during both steady-current and pulsed-current gas tungsten arc welding (GTAW). An analysis based on a wavelet transform was performed to extract the time-frequency content of the displacement signals obtained from numerical simulations. Our results confirm that the frequency of oscillations for a fully penetrated melt pool is lower than that of a partially penetrated melt pool with an abrupt change from partial to full penetration. We find that during transition from partial to full penetration state, two dominant frequencies coexist in the time-frequency spectrum. The results demonstrate that melt-pool oscillations profoundly depend on melt-pool shape and convection in the melt pool, which in turn is influenced by process parameters and material properties. The present numerical simulations reveal the unsteady evolution of melt pool oscillatory behaviour that are not predictable from published theoretical analyses. Additionally, using the proposed simulation-based approach, the need of triggering the melt-pool oscillations is expendable since even small surface displacements are detectable, which are not sensible to the current measurement devices employed in experiments.


Introduction
Molten metal behaviour during fusion-based welding and additive manufacturing affects energy transport in melt pools, which in turn influences their geometrical evolution [1,2] . Melt pool behaviour appears to substantially determine the properties, structure and quality of weldments or additively-manufactured products [3] . Correct control of melt pools ( e.g. penetration, shape and size) is therefore crucial to achieve desirable properties and to reduce the number of failures in production. However, melt pool control is usually a challenging task since its behaviour is characterised by highly non-linear responses to variations in process pa-material properties determine the fluid flow in a melt pool and can generate complex unsteady flow patterns as observed experimentally [7] and visualised numerically [8] . This fluid motion can result in self-excited oscillations of the melt-pool surface that may occur in various modes depending on melt-pool shape and size as well as the surface tension of the molten material [9][10][11] . Oscillations of the melt-pool surface can also be triggered by exposing the melt pool to external loadings such as welding current pulsation and filler droplet impingement. The melt-pool surface oscillations correlate with process conditions [12] and liquid metal properties [13,14] and can be utilised to measure, predict and control welding and additive manufacturing processes [15] .
The available literature on melt-pool surface oscillations, to a large extent, is experimental and focuses on employing the melt-pool oscillation signals for real-time penetration sensing [16] . The most common experimental techniques for detecting meltpool surface oscillations are based on the processing of arc voltage [15] , arc light-intensity [17] and laser vision [18] recorded during the process. These techniques only render information about the movement of the melt-pool surface and disregard the molten metal flow. Measurement schemes based on sensing the arc voltage or arc light-intensity cannot be employed to reflect local surface movements [19] since the recorded voltages and light intensities are averaged values across the melt pool at any given time during the process. Additionally, the signal-to-noise ratio in measurements decreases with a decrease in the amplitude of melt-pool surface deformations and may result in the obscuration of certain peaks in the frequency spectrum [20] . Numerical simulations, with an adequate level of accuracy, offer considerable potential to improve the current understanding of the pool surface oscillations by describing convection in the melt pool and its influence on geometrical evolution of the melt-pool boundary.
On the basis of a numerical study, Ko et al. [21] , 22 ] argued that the direction of the Marangoni flow affects the oscillations of the melt-pool surface during gas tungsten arc welding (GTAW). In their axisymmetric model, the surface tension of the molten metal was assumed to change linearly with temperature and the surfacetension temperature coefficient was assumed to remain constant ( i.e. ∂ γ /∂ T = C , where, C is a constant value). However, the surface tension of molten materials generally changes non-linearly with temperature if surfactants are present [23] , which can locally affect the direction of the Marangoni flow at the melt-pool surface and thus convection and subsequent solid-liquid phase transformations, resulting in complex three-dimensional flow patterns [3,8,24] . Cho et al. [25] implemented a three-dimensional model to investigate melt-pool oscillations during laser spot melting of titanium, assuming that the surface tension of molten titanium increases linearly with temperature, and stated that the period of oscillations changes during solidification. Huang et al. [26] employed a threedimensional model to study the melt-pool surface oscillations during pulsed-current GTAW and confirmed that the melt-pool oscillation frequency depends on the melt-pool size. Further investigations are required to foster our understanding of the influence of process parameters and material properties on melt pool oscillatory behaviour.
Deformations of the melt-pool surface can alter the surface power density and arc pressure distributions [21] and thus affect the melt-pool oscillatory behaviour. This effect is often neglected in numerical simulations implemented to study the meltpool dynamic behaviour, where the process is decoupled from the melt pool and is represented as a boundary condition at the gas-metal interface. Additionally, conventional Fourier transform analysis, which requires a sufficiently large sample size to return adequate stationary frequency spectra, is commonly employed to extract the frequency characteristics of melt-pool oscillations during fusion-based manufacturing processes [27] . This approach is not meticulously representative due to abrupt changes in the signals.
In the present study, a numerical approach is utilised to predict melt-pool formation, molten metal flow and oscillations of the melt-pool surface during stationary gas tungsten arc welding. Three-dimensional calculations are performed to study the effects of the process parameters as well as surfactants on meltpool surface oscillations. The model developed in the present work takes the influence of temporal surface deformations on arc powerdensity and arc pressure distributions into account, which is physically realistic and can affect the numerically predicted convection in the melt pool and hence the associated melt-pool surface oscillations. These effects are commonly neglected in published literature. Using the time-resolved data obtained from the numerical simulations, time-frequency analysis is conducted to enhance the understanding of melt-pool oscillatory behaviour. The present study elucidates molten metal oscillations during arc welding and arc-based additive manufacturing and offers an enhanced simulation-based approach to process development and optimisation. The proposed approach provides detailed information on molten metal flow in weld pools and the evolution of melt pool oscillatory behaviour.

Problem description
Gas tungsten arc welding of a stainless steel (AISI 304) plate with a thickness of H m = 3 mm , as shown schematically in Fig. 1 , is investigated numerically in the present work. Both steady current and pulsed current conditions (with different current pulsation frequencies) are studied. An electric arc plasma generates the power required to heat and melt the base material that is initially at an ambient temperature of T 0 = 300 K . The distance between the electrode tip and the nominal plate top-surface ( i.e. the arc length, ) is set to 2 . 5 mm before igniting the arc; however, arc length obviously changes as the melt-pool surface deforms. Changes in the arc length result in variations of the electric arc voltage U and thus the arc power and power-density distribution as well as the forces induced by the arc plasma [28,29] . These variations are accounted for in the present model by adjusting the source terms dynamically as described in Section 3 . The arc plasma is decoupled from the melt pool to reduce simulation complexity and run-time.
The computational domain in the present work is a rectangular cube that encompasses the base material and two layers of gas with a thickness of H a = 2 mm above and below the sample. The gas layers are included to monitor the melt-pool surface oscillations. The width of the computational domain in the x and y directions is W = 40 mm, which is considerably larger than the melt pool size. Variation of the surface tension with temperature is modelled using an empirical correlation that takes the influence of surfactants into consideration [23] . The thermophysical properties of the base-material and the gas employed in the numerical simulations are presented in Table 1 .

Mathematical model
A three-dimensional multiphase model is developed in the present work to predict the melt-pool dynamic behaviour and the associated heat and fluid flow. Both the molten metal and gas were treated as Newtonian and incompressible fluids. To describe heat and fluid flow in the melt pool and the associated solid-liquid phase transformation, the three-dimensional unsteady Navier-Stokes equations for conservation of mass, momentum and energy were cast as follow: where, ρ is density, u velocity vector, h sensible heat, t time, μ dynamic viscosity, k thermal conductivity, c p specific heat capacity at constant pressure, p pressure and H latent heat. Sum of the sensible heat h and the latent heat H returns the enthalpy of the material H and can be expressed mathematically as follows [32] : where, T is temperature, L latent heat of fusion, ψ local liquid volume-fraction and the subscript 'r' stand for reference condition.
The liquid volume-fraction ψ is assumed to correlate linearly with temperature [32] and was defined as follows: where, T s and T l are the solidus and liquidus temperatures, respectively.
The volume-of-fluid (VOF) method [33] was adopted in the present model to capture the position of melt-pool surface. In the VOF method, the scalar function φ represents the local volumefraction, and its advection was described the linear advection equation:

Dφ Dt
The value of φ ranges between 0 in the gas phase and 1 in the metal phase; hence the gas-metal interface locates in computational cells with 0 < φ < 1. Accordingly, the effective material properties in the model were determined based on the scalar function φ as follows: where, ξ corresponds to density ρ, specific heat capacity c p , thermal conductivity k and viscosity μ, and subscripts 'm' and 'g' stand for metal and gas, respectively.
The enthalpy-porosity technique [34] was employed to suppress fluid velocities in solid regions and to model the damping of fluid velocities in the mushy region, where solid-liquid phase change occurs between T s and T l . Accordingly, the sink term F d was added to the momentum equation and was defined as where, C is the permeability coefficient and a small constant number, equal to 10 −3 , to avoid division by zero. The value of the permeability coefficient C was set to 10 7 kg m −2 s −2 , determined from an earlier investigation on the sensitivity of numerical predictions to the value of the permeability coefficient [35] .
The continuum surface force (CSF) model [36] was employed to apply surface forces on the gas-metal interface, which considers surface forces as volumetric forces acting on the interface region. Hence, the source term F s was introduced into the momentum equation as follows: where, f s is the surface force per unit area, which includes surface tension, Marangoni and arc plasma forces and was defined as follows: where, γ is surface tension, κ the surface curvature κ = ∇ ·ˆ n , ˆ n the surface unit normal vector ˆ The correlation proposed by Sahoo et al. [23] was utilised to model the variations of the surface tension of the molten metal with temperature, which takes the effect of surfactants into account, and is expressed mathematically as follows: where, γ • m is the surface tension of pure molten-metal at the melting temperature T m , R the gas constant, s the adsorption at saturation, ψ an entropy factor, a s the activity of the solute, H • the standard heat of adsorption, and ( ∂ γ / ∂ T ) • the temperature coefficient of the surface tension of the pure molten-metal. Values of the properties used in Eq. (11) can be found in Sahoo et al. [23] . Variations of the surface tension of the molten metal with temperature for sulphur contents considered in the present work are shown in the supplementary materials.
In Eq. (10) f a is the arc plasma force, which comprises arc pressure f p (applied normal to the melt-pool surface) and arc plasma shear stress f τ (applied tangent to the melt-pool surface), The arc pressure f p was determined as follows [29] : where, μ 0 is the vacuum permeability equal to 4 π · 10 −7 H m −1 , R the radius in x -y plane ( i.e. , R = x 2 + y 2 ), I the welding current.
The distribution parameter σ p was determined based on the experimental results reported by Tsai and Eagar [28] for an argon arc with an electrode tip angle of 75 • as follows: where, is the local arc length in meters, I the current in ampere and σ p the distribution parameter in millimetre. Using this expression, variations of the arc pressure distribution due to changes in melt-pool surface shape were accounted for. F p was employed to negate changes in the total arc pressure force due to surface deformations through changes in ∇φ [37] , which may cause the total pressure force applied to the melt-pool surface ( fers from the available arc pressure force ( μ 0 I 2 /4 π ), and was defined as The term α was incorporated to match the theoretically determined arc pressure with experimentally measured values as suggested by Lin and Eagar [29] and Liu et al. [38] , and was defined as follows: The arc plasma shear stress f τ applies to the melt-pool surface in tangential direction, and was defined as follows [39] : where, the maximum arc shear stress τ max [40,41] , the arc shear stress distribution function g τ [42] and the surface unit tangent vector ˆ t [39] were defined as follows: Here, ¯ is the mean arc length over the melt-pool surface in meter, and r the position vector in the x -y plane. The distribution parameter σ τ was assumed to be a function of ¯ and welding current I and was approximated based on the data reported by Lee and Na [40] : The body forces applied to the melt pool F b include electromagnetic, thermal buoyancy and gravity forces. The model proposed by Tsao and Wu [43] was utilised to compute the electromagnetic force and the thermal buoyancy force was modelled using the Boussinesq approximation [44] . Accordingly, the three components of the total body force were defined as follows: where, g is the gravitational acceleration vector, β the thermal expansion coefficient, and σ e the distribution parameter for the electromagnetic force that is the same as σ p , according to Tsai and Eagar [28] . The heat generated by the arc was modelled through the incorporation of the source term S q in the energy equation and was defined as where, η is the process efficiency that changes linearly with welding current from 80% at 50 A to 70% at 300 A [45] . Variations of the arc voltage U with arc length and welding current was approximated as follows: where, U o is the electrode fall voltage equals to 8 V [46] , U e the electric field strength equals to 0.75 V m −1 [46] and U I the coefficient of variation of arc voltage with welding current equals to [47] . The distribution parameter σ q was determined using the data reported by Tsai and Eagar [28] as follows: To negate changes in the total heat input caused by melt-pool surface deformations, F q was defined as follows: Heat losses from the workpiece surface due to convection and radiation were accounted for by including S l in the energy equation and was determined as follows: where, K b is the Stefan-Boltzmann constant, E the radiation emissivity equals to 0.45 [48] and h c the heat transfer coefficient equals to 25 W m −2 K −1 [49] .

Numerical implementation
The model was developed within the framework of the proprietary solver based on the finite-volume approach, ANSYS Fluent [50] . The source terms in the momentum and energy equations as well as the temperature-dependent surface-tension model were implemented through user-defined functions. After analysing the sensitivity of numerical predictions to the computational grid size, the computational domain was discretised using non-uniform hexahedral cells with a minimum cell spacing of 35 μm in the melt pool central region close to the gas-metal interface. Cell spacing increases gradually towards the outer boundaries of the computational domain with a maximum cell spacing of 500 μm resulting in a grid with about 3 × 10 6 cells. The central-differencing scheme with second order accuracy was employed to discretise the diffusion and convection terms in the governing equations, and the transient advection terms were discretised using a first order implicit scheme. The PISO (Pressure-Implicit with Splitting of Operators) scheme [51] was utilised to couple velocity and pressure fields, and the PRESTO (PREssure STaggering Option) scheme [52] was employed for the pressure interpolation. The spatial discretisation of the gas-metal interface advection was performed using an explicit compressive VOF formulation [53] . A time-step size of 10 −5 s was chosen to achieve a Courant number ( Co = u t/ x ) less than 0.3. Scaled residuals of less than 10 −7 were defined as the convergence criterion for each time-step. Each case was run in parallel on 40 cores (Intel Xeon E5-2630 v4) of a highperformance computing cluster and took roughly about 500 h to complete. The accuracy of the present numerical simulations in predicting the evolution of the melt pool shape is validated against experimental data reported by Aucott et al. [3] , and the results are shown in Section 4 .

Time-frequency analysis
Melt-pool surface oscillations during fusion-based manufacturing processes are often time-variant due to the complex unsteady fluid flow in the melt pool and the associated melt-pool shape variations. Hence, the conventional Fourier transform (FT) analysis cannot reveal the non-stationary features of the signals whose frequency spectra varies with time [54] . To overcome the limitations of the Fourier transform and to access the time localisation of the frequency spectra, the continuous wavelet transform (CWT) was utilised [55] . The principle of the wavelet transform for signal processing is described comprehensively in Mallat [55] and is not repeated here. The Morlet wavelet, which is a Gaussian-windowed complex sinusoid, was employed as the mother wavelet that yields an adequate balance in both frequency and time domains and can adequately capture the coherence between harmonics. All timefrequency analysis in the present work was carried out in Python using the signals obtained from numerical simulations at a sampling rate of 10 5 Hz.

Validation of the numerical simulation
To validate the accuracy of the present model in predicting the melt pool shape during stationary gas tungsten arc welding, morphological evolution of the melt pool obtained from the present numerical simulation are compared to in situ experiments conducted by Aucott et al. [3] . In this problem, a gas tungsten arc welding with a steady-current of 125 A and an initial arc length of 1.027 mm was considered. The thermophysical properties of the alloy were obtained from Kidess et al. [8] . Fig. 2 shows the evolutions of the melt pool width and depth obtained from the present numerical simulation and the experimental data reported by Aucott et al. [3] , which indicates a reasonable agreement. The deviation between the numerically predicted melt pool shape and the reference data can be attributed to uncertainties associated with material properties employed in the simulation, particularly for temperatures above the melting temperature, and the simplifying assumptions made to develop the present model.

Results
In this section, the influence of welding current, current pulsation and surfactants on the melt-pool shape and its oscillatory behaviour during both partial and full penetration conditions are described. To monitor melt-pool surface oscillations, several monitoring points were defined in different locations over the melt-pool surface and displacement signals were recorded at these monitoring points. Signals received from a melt pool during fusion welding and additive manufacturing are generally time-variant and encompass abrupt changes. Therefore, time-frequency analysis based on a wavelet transform is preferred over the conventional Fourier transform analysis [55] . Signals recorded from the monitoring point located at the melt pool surface p c ( x, y, z ) = ( 0 , 0 , z surface ) in the period of 0.1 s to 3.7 s were utilised for time-frequency analysis and the results are shown in Fig. 3 . In case of steady-current GTAW with I = 110 A and sulphur content of 240 ppm ( Fig. 3 (a)), the melt-pool surface starts to oscillate with a fundamental frequency f 0 of about 195 Hz. Other frequencies are also observed in the time-frequency spectrum at f 1 ≈ 1.6 f 0 and f 2 ≈ 2.3 f 0 , which resemble the oscillation frequencies of a thin circular membrane predicted from theoretical analysis [20,56] . The frequency of oscillations decreases as the melt pool grows under partial penetration and reaches values of O(100) Hz. It is also found that the am-plitude of oscillations increases as the melt pool grows. Oscillations with a frequency of about 40 Hz then appear in the frequency spectrum at t ≈ 2.8 s, as indicated by arrow in Fig. 3 (a), when the melt-pool depth approaches the metal plate thickness ( i.e. full penetration). Interestingly, there is a period in the timefrequency spectrum ( 2 . 7 s < t < 3 . 5 s) that both high and low frequencies coexist while the transition from partial to full penetration is occurring. The oscillation frequencies obtained from the numerical simulations agree fairly well with the experimental measurements reported by Li et al. [19] and Yu et al. [27] and theoretical approximations reported by Xiao and den Ouden [12] , 15 ] for stationary GTAW. Fig. 3 (b) shows the displacement signals obtained from numerical simulations of GTAW with I = 150 A and sulphur content of 240 ppm and its corresponding time-frequency spectrum. The melting rate increases with increasing welding current from 110 A to 150 A, which is obvious because of the enhancement of total heat input to the bulk material and fluid velocities induced in the melt pool. Under partial penetration, the melt-pool width grows rapidly due to the outward Marangoni flow and the melt pool oscillates at a frequency of about 200 Hz up to t ≈ 1.75 s. The first arrow in the frequency spectrum of Fig. 3 (b) indicates that from t ≈ 1.8 s the frequency of oscillations decreases to values of about 125-140 Hz that persists for approximately 0.5 s, during which the increase in the melt-pool depth is more pronounced compared to its width. The second arrow in the frequency spectrum of Fig. 3 (b) at t = 2 . 2 s indicates that once the melt pool becomes fully-penetrated ( t ≈ 2.4 s), the melt pool oscillates at frequencies between about 70 Hz and 40 Hz. These observations demonstrate the significant effect of the melt-pool shape and size on oscillatory behaviour. Fig. 3 (c) shows the frequency spectrum for a pulsed-current GTA weld with a 240 ppm sulphur content. In this case, the melt pool size and shape, as well as flow pattern, changes continuously. During current peaks, the melt pool enlarges, and the melt pool surface depression increases due to augmented arc pressure. In the partial penetration state, after each current peak, the melt pool shrinks and the dominant frequency of oscillations increases with time and reaches values up to O(200) Hz. However, the dominant frequency of oscillations remains almost constant at about 60 Hz under a full penetration condition ( t > 2.8 s).
The behaviour of a 150 A steady current weld with a high sulphur (240 ppm) material ( Fig. 3 b) is compared with the same welding condition on a low sulphur (30 ppm) material ( Fig. 3 d). The melt pool evolves differently when the sulphur content of the material is lower. In contrast to the cases with 240 ppm sulphur content, after melt pool formation, the fluid flow on the melt-pool surface is directed outward forming a wide and shallow melt pool. Thermocapillary flow instabilities arose in this wide and shallow melt pool, generating an unsteady multi-cellular flow pattern in the melt pool [57] resulting in melt-pool surface deformations that are reflected in the displacement signals and the frequency spectrum, with frequencies of 200-310 Hz up to t ≈ 0.8 s. Changes in the melt-pool shape lead to the decay of these thermocapillary flow instabilities and the amplitude of melt-pool oscillations diminishes until t ≈ 2 s. Similar to the other cases described previously and as indicated by dashed lines in the frequency spectrum of Fig. 3 (d), the frequency of oscillations decreases as the melt pool grows under partial penetration. Flow instabilities grow again at t ≈ 2.3 s as the melt-pool depth increases, particularly in the central region, leading to an increase in the amplitude of surface oscillations. With further increase in the melt-pool depth, the frequency of oscillations drops to values of about 30-45 Hz at t ≈ 2.3 s, which is indicated by an arrow in the frequency spectrum of Fig. 3 (d).

Discussion
The results indicate that melt-pool oscillations strongly depend on melt-pool shape that oscillations are indeed affected by convection in the melt pool. Fig. 4 shows the melt-pool shape and temperature profiles over the melt-pool surfaces at different time instances during steady-current GTAW ( I = 110 A and sulphur content of 240 ppm). After arc ignition at t = 0 s, the bulk material heats up locally and a melt pool forms that grows over time. The forces applied to the molten metal in the melt pool ( i.e. Marangoni, electromagnetic, arc pressure, arc plasma shear and thermal buoyancy forces), drives the molten metal flow and gener-ates a complex three-dimensional flow pattern. This fluid motion transfers the heat absorbed from the heat source to the surrounding solid material. The Péclet number ( Pe = D u ρc p /k ), which quantifies the ratio of advective to diffusive heat transport, is greater than unity ( O(100) ) in the melt pool and indicates the remarkable influence of advection on total energy transfer and thus, on the melt-pool shape. Under partial penetration, the maximum temperature in the melt pool increases in time and reaches a quasi steady-state condition; however, it decreases slightly when the melt pool becomes fully penetrated, which in turn affects the average surface tension of the molten material. This decrease in maximum surface temperature is attributed to the flow of colder material from the bottom surface of the melt pool towards its top surface and the increase in the amount of heat loss from the surfaces due to radiation and convection.
Numerical predictions reveal an outward fluid flow from the central region of the melt pool that meets an inward flow from its outer edge. This occurs because of the change in the sign of surface-tension temperature coefficient ( ∂ γ / ∂ T ) at a specific temperature ( T ≈ 2115 K for an Fe-S alloy with a sulphur concentration of 240 ppm, according to the model proposed by Sahoo et al. [23] ) due to the presence of sulphur in the melt pool that is an active surface agent. The flow direction on the bottom surface of the melt pool is mostly inward since the temperature over the bottom surface is below the critical temperature of 2115 K. In this case, the maximum local fluid velocities predicted over the melt pool surface for partially penetrated melt pools are about 0.4 m s −1 , which agrees reasonably with the experimental measurements reported by Aucott et al. [3] . Nevertheless, the average bulk fluid velocity in fully penetrated melt pools are larger than those predicted in partially penetrated melt pools.
While the melt pool is partially penetrated, the melt-pool surface is depressed in the central region and is elevated in regions where the inward and outward streams meet each other. The unsteady interactions between these two streams coupled with continuous changes in the melt-pool shape disturb the fluid flow and consequently, the melt-pool surface oscillates with relatively small amplitudes and high frequencies. The amplitude of oscillations in-creases and the frequency of oscillations decreases as the melt pool grows. Once the melt pool becomes fully penetrated, the flow pattern changes due to Marangoni stresses induced over the lower melt-pool surface and changes in melt pool surface temperature. This change in flow pattern along with changes in the melt-pool shape and the average surface tension of the molten material result in an increase in the amplitude of melt-pool surface oscillations and a decrease in dominant oscillation frequency. The amplitudes of oscillations in a fully-penetrated melt pool appear to be larger in comparison with those of a partially penetrated pool. Changes in the melt-pool surface profile subsequently result in changes in the power available from the arc and its spatial distribution over the melt-pool surface, which enhances flow disturbances. Fig. 5 shows the evolution of the melt-pool shape and temperature distribution over its surfaces during steady-current GTAW with I = 150 A and sulphur content of 240 ppm. The average heat input to the bulk material increases with increasing welding current from 110 A to 150 A, resulting in a higher melting rate and a larger melt pool. The arc plasma force also increases with increasing the welding current, augmenting the melt-pool surface depression. Although the overall melt pool behaviour of the case with I = 150 A looks similar to that of the case with  ing welding with I = 150 A. Since the outer boundaries of the computational domain were assumed to be adiabatic, the melt pool keeps growing and collapses ( i.e. burn-through) eventually at t ≈ 3.8 s. Fig. 6 shows the evolution of melt-pool shape and temperature profile over its surfaces at different time instances during pulsed-current GTAW of AISI 304 with 240 ppm sulphur content. In this case, the current was pulsed from I b = 50 A to I p = 250 A at a frequency of 2.5 Hz and a duty cycle of 30% resulting in an average current of Ī = 110 A. During the peak current, the heat input and arc forces applied to the bulk material increase, resulting in an increase in the surface temperature. An outward fluid flow is observed at the melt pool top surface during peak current that makes the melt pool wider and enhances surface depression in the central region of the melt pool surface. The flow direction at the melt-pool surface reverses during the base current period due to the change in the sign of surface-tension temperature coefficient ( ∂ γ / ∂ T ) from negative to positive, which occurs because of melt-pool surface temperature reduction, resulting in a reduction in the melt-pool width and an increase in its depth. The melt-pool surface oscillates during the base current, and the amplitude of oscillations decays in time until the next current increase. In this paper, we have reported results for pulsed current welding at one specific pulsation frequency of 2.5 Hz. However, we have studied different current pulsation frequencies (1-10 Hz) as well, and the results are provided in the supplementary materials. All general qualitative behaviour discussed here were similarly observed for all studied pulsation frequencies.
The influence of Marangoni flow on melt-pool oscillations was examined by reducing sulphur concentration in the material from 240 ppm to 30 ppm, while the welding current I was set to 150 A (DC mode); the results are presented in Fig. 7 . In contrast to the cases with 240 ppm sulphur content, the fluid flow at the melt-pool surface is predominantly directed outward and transfers the heat absorbed from the heat source to the surrounding solid boundary. This leads to the formation of a wide and shallow melt pool, which differs notably from that predicted for the case with a sulphur concentration of 240 ppm (see Fig. 5 ). Additionally, this outward fluid flow enhances the melt pool depression in the central region. In steady-current GTAW, changes in powerdensity distribution over the melt pool surface due to surface deformations are more pronounced for the case with 30 ppm sulphur compared with those for the case with 240 ppm sulphur. This results in augmentation of flow disturbances because of changes in temperature distribution and thus Marangoni stresses induced over the melt-pool surface. These effects lead to oscillations with amplitudes that are larger than those predicted for the case with 240 ppm sulphur content.
The melt-pool oscillation frequencies obtained from the simulations are compared with the experimental data reported by Xiao and den Ouden [15] and theoretical predictions calculated using the model proposed by the same authors [12,15] for both partial and full penetration conditions, and the results are shown in Fig. 8 . Xiao and den Ouden [15] employed a plate made of AISI 304 with 170 ppm sulphur content and they applied a pulsed current with a base current I b ranging between 70 A and 150 A and a peak current I p ranging between 100 A and 180 A at a frequency of 25 Hz under partial penetration and 10 Hz under full penetration. The first two modal frequencies of oscillations under partial penetration are calculated from the following equations [15] : where, D e is the equivalent melt-pool diameter and equals to the melt-pool diameter on the top surface D t . The frequency of oscillation under full penetration is calculated from the following equation [15] : where, D e is the equivalent melt-pool diameter under full penetration condition defined as follows [15] : Here, D t and D b are the melt-pool diameter measured on the top and bottom surfaces, respectively. To calculate oscillation frequencies from the theoretical models, the average surface tension of the molten metal γ is assumed to be 1.6 N m −1 , which is the mean value of surface tension in the temperature range of 1673-2300 K for an Fe-S alloy with a sulphur concentration of 240 ppm, according to the model proposed by Sahoo et al. [23] .
The results indicate that depending on the processing condition, the melt pool can oscillate in different modes. It seems that the melt pool is more likely to oscillate in the sloshy mode (mode 2) with a reduced mean welding current. The predictions obtained from the theoretical models agree reasonably with the present numerical results and the experimental measurements, which confirms their adequacy in predicting the melt pool oscillation frequencies. The deviation between theoretical predictions and the numerical and experimental data is attributed to the simplifying assumptions made to develop the models and uncertainties in prescribing the average surface tension of the molten metal. It should be noted that although theoretical models can provide information about oscillation frequencies, they fail to predict the variation of oscillation amplitude as well as changes in the oscillation mode, which are influenced by the complex melt flow. Additionally, we have shown that by taking the proposed simulation-based approach the need of triggering the melt-pool oscillations by a current pulsation, mainly to amplify the oscillation  amplitudes to make them sensible for measurement devices, is expendable.

Conclusions
A simulation-based approach was implemented to study heat and fluid flow in molten metal melt pools and associated surface oscillations during gas tungsten arc welding (GTAW). Utilising the proposed approach, the effects of surfactants and various welding process parameters on melt-pool oscillations were investigated. Time-frequency analysis based on a wavelet transform was carried out to improve our understanding of melt pool oscillatory behaviour with time dependent responses that may encompass abrupt changes. The results obtained from the present numerical simulations are compared with experimental and theoretical results. Using this approach, the frequency and amplitude of melt pool oscillation and changes in the oscillation modes were predicted, which are not predictable using published theoretical models and are generally difficult to measure experimentally.
The melt-pool oscillations strongly relate to melt-pool shape that is determined by convection in the melt pool to a large extent. Convection in the melt pool and thus the melt-pool shape are both affected by flow instabilities that arise rapidly in the melt pool. The spatial distribution of arc power-density and forces imposed on the melt pool change with surface deformations that can affect the temperature distribution over the melt-pool surface, and thus flow instabilities that are often dominated by Marangoni flow. This effect, which is often neglected in numerical simulations of melt pool behaviour, is accounted for in the present study. Oscillations with relatively low amplitudes and high frequencies were observed during partial penetration. The frequency of melt-pool oscillations decreases as the melt pool grows. Our numerical predictions confirm the existence of a sharp drop in oscillation frequency when the melt pool becomes fully penetrated. We find that there is a period that both high and low frequency coexist while the transition from partial to full penetration is occurring. Additionally, the results show that the melt pool oscillation mode can change during the process depending on material properties and prescribed process parameters. The melt pool oscillation and thus the process stability can be optimised through adjusting welding process parameters and adjusting the concentration of surfactants.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.