Numerical study on the thermodynamic behavior of de-icing liquid droplets impacting walls

Spraying de-icing fluid is a key method to ensure the safe operation of aircraft in icy and snowy weather. The film aggregation and internal mixing of de-icing fluid droplets on the aircraft skin during a collision are crucial. Considering the rheological properties of the molecular viscosity change of the de-icing fluid droplets during the collision and the heat transfer model of the heat loss after the impact, the phase field method is used to capture the gas–liquid interface, and a thermal pressure/viscous coupling model is constructed. The thermodynamic behavior of different axial distances is calculated. The results show that, as the dimensionless axial distance of the droplet increases, the spreading length of the fused droplet decreases instead, and the heat transfer rate of the droplet increases with the increase in spreading length. After stabilizing, the increase or decrease in the heat transfer rate depends on the strength of the heat transfer between the liquid layers. As the dimensionless axial distance increases, the internal flow of the droplet weakens and, between the droplet and the wall, the heat flux density gradually decreases and the average temperature drop of the droplet becomes gradual.


Introduction
Ice and snow constitute substantial threats to the safe operation of aircraft during winter.To mitigate these challenges, prevalent de-icing approaches encompass electrical heating de-icing [1][2][3] and the application of de-icing agents [4].Jiang and Wang [1] scrutinized the effectiveness of the electric pulse de-icing method through a combination of numerical simulations and experimental investigations, unveiling superior performance achieved with a dual-coil de-icing system.Mohseni and Amirfazli [2] introduced an embedded heating element tailored for wings made of composite materials, providing empirical evidence of its viability.However, the application of electrical heating de-icing is predominantly tailored for specific components such as wings, prompting the practical integration of de-icing liquids in real-world scenarios.
Type II, III, and IV de-icing fluid has garnered attention for its shear-thinning non-Newtonian properties that offer excellent anti-splash and anti-freezing performance [5].These attributes render non-Newtonian de-icing fluids more effective in adhering to the aircraft skin surface, leading to superior de-icing and anti-icing outcomes.The utilization of non-Newtonian de-icing fluid spray involves a myriad of physical phenomena, encompassing droplet merging [6][7][8][9], fragmentation [10], and turbulent flow [11].An in-depth investigation into these phenomena contributes to a comprehensive understanding of the efficacy of non-Newtonian de-icing fluid on aircraft skin.Furthermore, such studies facilitate the enhancement of its efficiency, ensuring heightened levels of traffic safety.
Substantial advancements have been achieved in the exploration of non-Newtonian droplet collision characteristics.Experimental investigations constitute a pivotal research approach for comprehending the intricacies of shear-thinning droplet collisions.An and Lee [12,13] conducted insightful experiments focusing on the behavior of xanthan gum solution droplets impacting solid walls.Their findings underscored the profound influence of shear-thinning characteristics on droplet dynamics, notably accelerating the spreading stage while concurrently suppressing the retraction process.German and Bertola [14] delved into the impact behavior of xanthan gum droplets on hydrophobic walls, suggesting that, when describing shear-thinning rheological properties using the power law model, the shear-thinning index exerts minimal influence on the droplet spreading process, with the dominance of the consistency coefficient.In a study by Laan et al. [15], blood droplets characterized by shear-thinning properties were investigated, revealing that the high shear rate viscosity solely affects the maximum spreading diameter of the droplets.
Concurrently, numerical simulation methodologies have become widely adopted in the exploration of shear-thinning droplet collisions and de-icing phenomena.Hansman et al. [16] conducted observations on the formation of feather-like ice in an icing tunnel, aligning with the dry haze accretion model proposed by Personne.Ferro et al. [3] made substantial strides by introducing a comprehensive CFD icing model, thereby contributing significantly to the realm of de-icing assessment.Bragg et al. [17] further delved into the repercussions of ice accretion on the aerodynamic parameters of aircraft, suggesting a methodology for detecting icing through the evaluation of changes in steady-state aircraft parameters.In a study by Xia et al. [18], macro-scale multi-phase flow CFD simulations were employed to quantitatively analyze the fragmentation threshold of droplets impacting rough surfaces.Sun et al. [19] utilized the meso-scopic lattice Boltzmann method to investigate the merging and internal mixing of two non-Newtonian droplets colliding with each other.Nomura et al. [20] applied the micro-scopic MPS (moving particle semi-implicit) method to scrutinize the rupture of droplets enveloped by a gas film.Li et al. [21] proposed a novel definition of the Ohnesorge number, incorporating the characteristic shear rate induced by air flow on droplets.Their findings indicated that shear-thinning properties alter the local apparent viscosity of droplets, thereby influencing fragmentation morphology.
Considerable achievements have been made in the study of the impact of individual Newtonian liquid droplets on solid surfaces [22][23][24].The research on non-Newtonian double droplets, however, has primarily focused on the airborne collision stage [19,25], while this paper places more emphasis on investigating the behavior of droplets colliding and merging on solid surfaces.The viscosity of non-Newtonian de-icing fluids fluctuates with the shear rate, introducing complexity to the droplet behavior during impact.Numerical simulation technology serves as a valuable tool for computationally replicating the droplet collision process, thereby enhancing comprehension of the efficacy of non-Newtonian de-icing fluids in practical scenarios.This technology not only facilitates a deeper understanding of the non-Newtonian de-icing fluid behavior but also offers essential insights for its practical application.Consequently, numerical research on the collision dynamics of non-Newtonian de-icing fluid droplets holds the potential to furnish robust support for its application, contributing to the secure and sustainable advancement of air traffic.

Mathematical model
The challenge of high-temperature de-icing fluid droplets impacting aircraft skin can be conceptualized as the interaction between individual droplets and the aircraft skin.While numerous studies have explored the dynamics of single droplet and wall collisions, this paper concentrates on the intricacies of double-droplet collision on a wall.The Euler-Euler multi-phase flow model, employing the volume of fluid (VOF) method along with the energy equation and the Navier-Stokes equation, is employed to investigate the physical process of film fusion and internal mixing of non-Newtonian droplets during collisions.Additionally, the accompanying behaviors of energy, momentum, and viscosity loss are scrutinized.In this section, the process is streamlined, and a two-dimensional droplet collision simulation model with a solid wall is implemented using Fluent software.The particulars of this process are outlined as follows: a circular droplet with a diameter D 0 , initial velocity U 0 , and initial temperature T 0 collides with a droplet of equal volume attached to the wall under air and gravity conditions.The distance between the axes of the two droplets is denoted as " l " and the vertical distance is "h = 1.5 mm." The adhered droplet is formed under the same conditions as the droplet/wall collision scenario, serving as the initial condition for subsequent calculations.The length of the adhered droplet is 2 mm, and its thickness is 0.4 mm.The wall material is aluminum, with a thermal conductivity of 202.4 W/(m•K), initial temperature T wall and thickness h wall .As the model in this article is two-dimensional, the conclusion of eccentric collision only applies to two colliding liquid cylinders, not droplets.The heat transfer between the droplet, the wall, and the air is considered, and the impact process's dynamic characteristics and changes in physical parameters are examined, with a detailed analysis of the influencing factors.The numerical model is depicted in Fig. 1.

Governing equations
Assuming that the multi-phase flow is incompressible and immiscible, the Navier-Stokes equations for mass and momentum conservation are as follows: (1) where U (U = (u, v, w)) is the velocity vector, t is the time, ρ is the fluid density, and g is the gravitational acceleration.The VOF model is used to track the gas-liquid interface, and the tracking of the interface between phases is accomplished by solving a continuity equation for the volume fraction of one (or more) phase(s).For the q th phase, this equa- tion has the following form: where ṁpq is the mass transfer from phase p to phase q , ṁqp is the mass transfer from phase q to phase p , U q is the velocity of fluid q , ρ q is the density of fluid q , and S q α is the mass source.The volume fraction of the primary phase will be calculated according to the following constraint: where α q is the volume fraction of phase q , and there are three possible cases: Therefore, the density ρ and the viscosity µ of the entire computational domain are con- sidered as an effective fluid with physical properties, and are obtained by the following method: (3)    α q = 0; The cell is empty (of the fluid).0 < α q < 1; The cell contains the interface between the fluid and one or more other fluids.α q = 1; The cell is full (of the fluid).
( where the subscripts L and G refer to liquid and gas, respectively.Since the heat transfer between the droplet and air and the wall is involved, the energy model is also added: where E is the energy, T is the temperature, S h is the source term, k eff is the effective thermal conductivity ( k + k t , where k t is the turbulent thermal conductivity defined by the turbulence model used), j j is the diffusion flux of substance j , h j,q is the enthalpy of substance j in phase q , and j j,q is the diffusion flux of substance j in phase q .The first three terms on the right side of the equation represent the energy transfer due to conduction, species diffusion, and viscous dissipation, respectively.The VOF model treats energy as a mass-averaged variable: Surface energy is a physical property that represents the energy required to form a new unit area.The surface tends to reduce its own area to minimize energy, which is the reason for the formation of spherical droplets and bubbles, as well as the surface tension of liquids.The formula for the surface energy can be expressed as: where γ is the surface tension and A is the surface area.

Boundary condition
The initial conditions for the two-dimensional phase field model are illustrated in Fig. 1.
The dimensions of the computational domain are determined by the spreading diameter and rebound height of the impacting droplet.A droplet with a diameter of D 0 is positioned at a specific distance from the bottom of the solid surface.If the impact of the droplet transpires in an unbounded domain, the upper and lateral boundaries of the three-dimensional domain are designated as open flow boundaries.The initial velocity of the colliding droplet is imposed through a volume force.In accordance with the theory postulated by Kistler [26], it is hypothesized that the dynamic contact angle ( θ d ) of the moving contact line is a function of the three-phase contact line velocity: where f H and f −1 H are the Hoffman function and its inverse function, respectively, whose functional descriptions are as follows: where x is the wetting radius of the impacting droplet and C a is the capillary number, which is described as follows: where σ is the surface tension and V CL is the moving contact line velocity, which can be equal to the velocity at the interface of the first computational point above the wet wall, which is described as follows: when V CL > 0 , the droplet is in the spreading phase.In the receding phase, its value is negative ( V CL < 0 ).θ e is the equilibrium contact angle, which is set to 40°.

Material properties
The simulated de-icing fluid droplet exhibits shear-thinning characteristics, with its viscosity diminishing as the shear rate increases.The rheological characteristics of the droplet can be aptly described by the power law model [5]: where γ is the shear strain rate, K is the consistency coefficient, and m is the power law index.Considering the droplet impact conditions, the characteristic shear rate of the droplet (U/D 0 ) is of the order of 10 3 s −1 .To avoid the distortion of the power law model when the shear rate is too small, the power law model is truncated and corrected, and the viscosity when the shear rate is less than 25s −1 is regarded as constant, which is equal to the zero-shear viscosity ( η 0 = 3 kg(m • s) ).Since the viscosity of the shear-thin- ning fluid is a variable value, the generalized Reynolds number Re n is used to character- ize the state of the shear-thinning droplet impacting the wall.The general expression for the Reynolds number of Newtonian fluids is: where ρ is the fluid density, D 0 is the characteristic diameter, U is the fluid velocity, and µ is the dynamic viscosity of Newtonian fluids.Use a power-law model to describe the viscosity of shear-thinning fluids, which is: where K is the dynamic viscosity, γ is the shear rate, γ0 is the reference shear rate, and m is a non Newtonian exponent.Substitute the non Newtonian viscosity model into the Reynolds number expression of Newtonian fluids: algebraically simplify the above equation to obtain the generalized Reynolds number: The other parameters of the de-icing fluid droplets are shown in Table 1.Given the impact of natural convection in air and the corresponding density variation with temperature, the Boussinesq model is adopted for density, with the remaining parameters set to their default values.As the aircraft skin predominantly comprises aluminum alloy, the wall material is appropriately designated as aluminum to align with real-world conditions.

Model validation
The numerical simulation has previously demonstrated successful validation for the collision and coalescence of Newtonian fluids [24].To further substantiate the numerical methodology, validation is extended to encompass shear-thinning non-Newtonian droplets.The validation was executed using shear-thinning Type II deicing fluid with a contact angle of 40°, maintaining consistent initial parameters outlined in Table 1.The viscosity was characterized using the power law model: η( γ ) = K γ m−1 , with the consistency coefficient K set to 2.2 and the power law index n set to 0.42 (Fig. 2b).Comparative analysis between experimental and simulation outcomes was conducted across various dimensionless times ( τ = tU 0 /D 0 ).The two sets of results exhibited fundamental consistency in morphology during the spreading phase (Fig. 2a).Furthermore, a comparative assessment was undertaken regarding the temporal variations in spreading diameter and central liquid film thickness during the droplet fusion process (Fig. 2c).The simulation error of the spreading diameter is 12%, and the simulation error of the central liquid film thickness is 6%, falling within (20) Table 1 The physical properties of de-icing fluid droplets [5] Parameters The accuracy and convergence of the model are closely related to the mesh delineation.While the bubbles generated at the liquid-liquid interface during the droplet/liquid droplet collision are observed, to capture this phenomenon a fine mesh is required.This study employs fine grids and high-order numerical schemes to accurately capture interface dynamics and liquid film deformations.An adaptive grid is implemented to refine the gas-liquid phase interface, as illustrated in Fig. 3. Adaptive grid: Refine when the proportion of liquid in the cell is greater than 0.5, coarsen when the proportion is less than 0.1, and generate a grid every ten steps.The larger refinement level is 1, the number of encryption layers is 3, and the saved data type is double precision.Additionally, a gradient surface tension model is incorporated to simulate the interactions between the gas and liquid phases, including effects such as the Marangoni effect induced by temperature variations and the expansion of liquid bridges during the mutual contact of two droplets [26].

Results
The double-droplet collision is affected by a series of impact parameters, such as surface characteristics (roughness, contact angle θ between gas and liquid phases), droplet physical parameters (droplet density ρ l , droplet viscosity η , surface tension σ between droplet and gas phase, etc.), and impact parameters (initial droplet diameter Fig. 2 The validation of the Type II de-icing fluid model includes a a comparison between experimental and simulated results for the collision process, b power-law fitting of the rheological model, and c a qualitative comparison between experimental and simulated results for the impact process D 0 , initial impact velocity U 0 , etc.).In addition, a series of dimensionless impact parameters are customized, as shown in Table 2.

The effect of shear thinning
The initial investigation focused on the collision of two droplets with shear-thinning properties on a wall, where one droplet was attached to the wall and the other droplet had a velocity U 0 with a dimensionless axial distance B equal to 0. Figure 4 illustrates the contour variation, velocity field, shear rate (left), and laminar viscosity (right) of the double droplets during the collision process with time τ .The maximum generalized Reynolds number ( Re n ) for the colliding droplets was 0.64, significantly below 100, indicating that viscous forces on spreading could not be ignored [27][28][29].Moreover, at the initial collision stage, the Ohnesorge number ( Oh ) was calculated as 2.13, signifying that inertia plays a more substantial role in droplet motion relative to viscosity and surface tension forces as Oh > 1 .As the collision progresses, the droplet's viscosity gradually increases, resulting in a corresponding rise in the Ohnesorge number.This observation suggests that, throughout the entire motion process, the droplet's shape and fusion are more influenced by inertial forces.At τ = 0.4 , the weaker resistance of the attached droplet to the wall led to the deformation of the merged part, while the upper half remained spherical.The shear rate cloud map indicated that the viscosity-affected droplet produced shear at the liquid-liquid interface due to resistance.The shear rate of the merged droplet increased gradually from the center axis to the outside, reaching a maximum shear rate of γ = 10548s −1 .The lam- inar viscosity of the entire droplet remained low at this point, close to the infinite shear viscosity ( η ∞ = 0.5 Pa • s ), still orders of magnitude larger than water viscosity, resulting in significant energy dissipation due to viscosity during the collision process.Additionally, it was evident that the shear rate at the liquid film surface surpassed the shear rate within the bulk, potentially triggering the formation of capillary waves and vortices near the liquid film.These vortices and capillary waves could contribute to increased shear rates near the surface.At τ = 0.8 , under the impact of the Fig. 3 Adaptive encryption mesh colliding droplet, a velocity field was generated inside the spreading droplet, gradually increasing along the liquid-liquid contact line to both sides, resembling a stratified flow.Over time, at τ = 1.2 − 1.6 , the shear rate outside the merged droplet began to decrease but remained at a high order of magnitude around 10 s −1 .The recovery of laminar viscosity was not pronounced until τ = 6.8 , when the laminar viscosity of the merged droplet gradually recovered from the bottom and edge portions.Concurrently, during the liquid film fusion process, the viscosity around the bubbles between liquid layers significantly decreased compared with other regions.The presence of bubbles appeared to enhance the local liquid's fluidity.Since bubbles represent a discontinuous phase in the fluid, their existence may mitigate overall internal frictional effects, facilitating local liquid flow and ultimately increasing the local shear rate, resulting in lower viscosity around the bubbles.
The viscous force acting on the droplet surface significantly influences the deformation and coalescence process when droplets collide.Therefore, the time variation of the average laminar viscosity ( η ) of the Newtonian and non Newtonian double droplets during the collision process was investigated (Fig. 5).It can be seen that the viscosity of Newtonian droplets remains at 3 Pa•s during collision.Higher viscosity greatly reduces the fusion rate of droplets and also suppresses their spreading on the wall.However, when the two shear thinning droplets come into contact, η is at a relatively low level, approxi- mately equal to its infinite shear viscosity ( η ∞ ).As the droplets spread and coalesce, part of the kinetic energy of the impacting droplet is converted into the surface energy of the coalesced droplet, and another part is dissipated by viscous forces and surface friction forces.As a result, the droplet motion decelerates, the shear rate decreases, and the laminar viscosity gradually recovers from the contact area and the edge with the wall.Eventually, it can fully recover to zero shear viscosity ( η 0 ).However, in practice, de-icing fluid may experience viscosity loss due to factors such as a decrease in the molecular weight of the thickening agent resulting from prolonged shear and high temperature exposure [5].This paper does not consider the impact of physical and chemical factors on molecular weight.If a more in-depth analysis of viscosity loss is required, the numerical model would need further modification.
It is important to note that a series of bubbles was generated during the collision process of the high-viscosity shear-thinning double droplets (Fig. 6).This is consistent with the bubble trapping phenomenon in methanol aqueous solutions studied by Zhang et al.

Dimensionless Parameters Expression Descriptions
Weber number We = ρ l U 0 2 D 0 /σ Measuring the relative importance of inertial force and surface tension The ratio of the axial distance between two droplets to the initial diameter [30].The compression of gas inside the high-speed impacting droplet may exceed the surface tension that encapsulates it, leading to outward bursting and the formation of bubbles ( τ = 0.4 ).Simultaneously, the high-laminar viscosity of the droplet hinders

The effect of the dimensionless axial distance
The coalescence and spreading of the droplets under different dimensionless axial distances ( B = 0 − 1 ) were investigated, with a primary focus on the viscosity and mor- phology changes of the two droplets during the coalescence process (Fig. 7).Contrary to a simple increase in the dimensionless axial distance B of the two droplets, the final spreading length after coalescence gradually decreases.When the two droplets come into contact, the inertial force and shear stress between the moving droplets induce fluid flow in the contact area, reducing the viscosity of the droplet contact area ( τ = 0 ) and lessening the resistance to droplet motion.Propelled by inertial force, the droplets move to both sides.The shear stress typically transmits from the contact point between the two droplets and gradually weakens towards both sides.The viscosity reduction ( τ = 0 ) caused by shear effect also weakens accordingly.As the dimensionless axial distance B increases, the contact area between the two droplets decreases, and the transmission of shear stress in the contact area weakens.The propagation of the inertial force along the spreading direction gradually diminishes, and the fluid at a greater distance from the impact center area experiences greater resistance, ultimately resulting in a decrease in the spreading length of the coalesced droplet.The variation of the line contact points ( X L and X R ) and the center ( X C = (X L + X R )/2 ) along the X direction of the coalesced droplet after coalescence was further analyzed Fig. 6 The variation of the pressure field during the collision process (left: overall view, right: local view) (Fig. 8).The change in X R is larger than that in X L , being attributed to the higher shear stress gradient at the right impact point during the impact.The viscosity near the right contact line decreases, reducing the resistance to motion and leading to a significantly faster movement rate of the right line contact point than the left one.Simultaneously, as the dimensionless distance B increases, the growth rate of X R accelerates, while the growth rate of X L slows down, eventually reaching a stable state.Despite differing rates of change, the stable value of X C consistently hovers around X = 0.13 , which can be approximately expressed as X R = X L + l + a .Here, X R is the displacement of the right (close to the impact point) line contact point, X L is the displacement of the left (far from the impact point) line contact point, l is the axial distance between the two Fig. 7 The time sequence phase diagrams of two droplets' collision at different distances (black arrows indicate the velocity direction) Fig. 8 The variation of the line contact point along the X direction of the coalesced droplet with time droplets, and a is a constant.In this study, a is determined as 0.13 for the shear-thinning droplets under investigation.

The heat transfer during impact process
The de-icing fluid is usually heated to a specified temperature during use, and its viscosity does not strictly follow the rule of decreasing with increasing temperature [5].Therefore, the heat transfer characteristics of the binary droplet collision process will be further studied.
Neglecting the heat loss due to droplet phase change and evaporation, the heat transfer during high-temperature binary droplet collision with a wall mainly includes thermal conduction between the droplet wetting area and the wall and convective heat transfer between the droplet and the surrounding air.The temperature field of the binary droplet collision process is shown in Fig. 9.When the two droplets contact ( τ = 0 ), heat transfer has already started, and the temperature gradient near the droplet/wall contact line is large ( 1308 K/mm ).As the two droplets merge and spread further ( τ = 0.4 ), the tem- perature at the droplet edge begins to decrease, and the temperature of the contact part between the droplet and the wall drops to 341 K , while the temperature of the central part of the droplet ( 352 K ) does not change much.Due to the viscous effect between the droplet and the air, the air "slip layer" on the droplet surface drives the gas to diffuse to both sides ( τ = 0.8 ), spreading the heat out.When τ = 4.8 , the droplet spreading reaches its maximum diameter, and the heat exchange between the droplet and the wall gradually increases.The temperature gradient in the center of the wetting part of the droplet is large ( 219 K/mm ), and the maximum temperature gradient of the contact part between the liquid film and air is located at the center of the upper surface ( 728 K/mm ).As time increases ( τ = 14.4 ), the air on the upper surface of the liquid film is heated and collides, reducing its density, and natural upward convection takes away heat.
The shear-thinning property reduces the flow resistance of the collision part, and the internal flow velocity of the droplet affects the heat transfer intensity of the droplet.The heat transfer of the droplet with different dimensionless axial distances B was studied.The temperature cloud map and streamline cloud map of the binary droplet/wall collision process are shown in Fig. 10.As the dimensionless axial distance B increases, the droplet spreading length gradually decreases, and the rate of temperature decrease of Fig. 9 Temperature field variation during the collision process the droplet also becomes slower.Due to the heat transfer between the wall and the air, the temperature gradient generated by the upward movement inside the droplet, driven by the surface tension and hindered by the high viscosity area of the droplet (Fig. 4, τ = 6.8 ), eventually leads to the internal streamline of the droplet rising from the middle to both sides.A significant vortex is generated in the low-viscosity area inside the droplet, and the increase in flow velocity leads to enhanced heat transfer inside the droplet.At the same time, the proportional reduction in the spreading length of the droplet after collision also directly leads to a gradual decrease in the rate of temperature decline of the droplet.
The average droplet temperature and wall heat flux density with different dimensionless axial distances B were analyzed (Fig. 11).It was found that, as time increased, the heat flux density q between the droplet and the wall increased with the spread- ing length until stabilization occurred when the heat transfer between the liquid layers was equal to that between the droplet and the wall, while the average droplet temperature maintained a uniform decreasing speed.As the droplet axial distance B increased, the heat transfer from the droplet to the wall was reduced.The aver- age temperature gradually rises, and the rate of decrease becomes slower.At B = 1, the droplet experiences the smallest average temperature drop, reaching its maximum value, and the corresponding heat flux density showed the opposite trend.The decrease in the average droplet temperature was mainly affected by natural air convection and wall heat transfer.During the spreading stage, as the contact area with the wall decreased, the heat flux density between the droplet and the wall rapidly decreased, and the rate of temperature decrease of the droplet also becomes slower.When the spreading length reached its maximum, the droplet heat flux density began to decrease, and the average droplet temperature maintained a uniform decreasing speed.At the same time, due to the overall temperature decrease in the droplet, the internal temperature gradient was small, which weakened the internal convection and eventually led to a weakening of heat transfer with the wall.Changing the axial distance B of the droplet collision will significantly affect the heat transfer characteristics of the droplet.The temperature loss at B = 1 was minimal, about 2.7 times that at B = 0.

Conclusions
The phase field method was employed to track the gas-liquid interface, and the power law model was integrated to describe the shear-thinning characteristics of the droplets.The study considered heat transfer during droplet collision processes and conducted a numerical investigation into the collision and fusion of two droplets on a wall surface.However, the use of a two-dimensional model neglected the true threedimensional motion trajectory of the droplets.Additionally, the size effects inherent to CFD prevented the micro-scopic representation of the liquid film coalescence process.Future work may explore employing a three-dimensional molecular dynamics approach to further simulate the micro-scopic scale interactions of forces.The current research yields the following conclusions: (1) The collision and coalescence of Type-II de-icing fluid droplets exhibit complex dynamic effects due to the high viscosity at zero shear rate and the local viscosity reduction caused by the shear rate increase.The bubbles trapped during the collision form a liquid-liquid interface, resulting in a stratified flow of the droplet, which impedes the internal mixing of the droplet.Meanwhile, the high viscous dissipation suppresses the droplet retraction.By considering the viscosity, flow characteristics, and thermal conductivity properties of the droplets, the formulation and performance of de-icing fluids can be optimized, thereby improving the de-icing efficiency of the de-icing fluid.(2) As the dimensionless axial distance B of the droplet center increases, it is observed that the contact line velocity on the right side (near the impact point) is significantly faster than that on the other side.The contact line displacements on both sides satisfy the relation X R = X L + l + a , indicating that the displacement increment of the right contact line is equal to the displacement increment of the left contact line plus the initial axial distance between the two droplets and a constant k.The final spreading length ( X R − X L ) decreases with the increase in the dimensionless axial distance B.
(3) This study takes into account the heat conduction between the droplet and the wall during the collision process, along with the natural convection of air.With increasing time, the heat flux density q between the droplet and the wall increases with the spreading length, stabilizing when the heat transfer between the liquid layers and the heat transfer between the droplet and the wall are in equilibrium.Concurrently, the average temperature of the droplet undergoes a uniform decrease.As the dimensionless axial distance increases, the heat flux density from the droplet to the wall decreases proportionally, and the average temperature drop becomes gradual.
Studying the heat transfer process contributes to understanding energy consumption and heat losses.By optimizing the heat transfer performance of de-icing fluids, energy consumption can be reduced, and carbon emissions can be lowered, thus making it more environmentally friendly.
Re n = ρ l D 0 U 0 (2−m) /K Measuring the relative importance of inertial and viscous forces Ohnesorge number Oh = η l / √ ρ l σ D 0 Measuring the relative importance of the relationship of surface tension, viscosity, and inertia Dimensionless time τ = tU 0 /D 0 The ratio of the product of time and speed to the initial diameter Dimensionless axial distance B = l/D 0

Fig. 4 Fig. 5
Fig. 4 Time series cloud map of collision between Newtonian and non Newtonian droplets.a is the cloud map of a non Newtonian droplet collision, where the left is the shear rate cloud map, and the right is the viscosity cloud map; b is the volume fraction map of Newtonian droplets, and the black arrow represents the velocity vector

Fig. 10
Fig. 10 Temperature field and streamlined cloud map of two droplets' collision with different dimensionless axial distance

Fig. 11
Fig. 11 Temperature and heat flux density variation during the collision process (a average droplet temperature, b wall heat flux density)