Numerical Model Study of Multiple Dendrite Motion Behavior in Melt Based on LBM-CA Method

: In this paper, a new method is proposed to solve the solute ﬁeld of moving grains, and a Cellular automaton (CA)-Lattice Boltzmann method (LBM)-Semi rebound format(Ladd) coupling model which can accurately simulate the motion behavior of multiple dendrites is established. The growth process of microstructure in the solidiﬁcation process of Al-4.7% Cu alloy ingot was calculated by Cellular automaton (CA) method, the momentum, heat, and mass transfer processes were calculated by Lattice Boltzmann method (LBM), and the melt-dendrite sharp interface interaction was treated by Ladd method. The reliability of the model is veriﬁed, and then the growth and movement of single dendrite and multiple dendrites under the action of gravity ﬁeld are simulated. The simulation results show that the growth and movement mode of multiple dendrites are quite di ﬀ erent from that of single dendrite, which is shown in two aspects: (1) the original motion state of dendrites is changed by the combination of ﬂow ﬁeld, which slows down the falling speed of dendrites to a certain extent; (2) the fusion of solute ﬁeld between dendrites changed the original growth mode of boundary dendrites and increased their rotation speed.


Introduction
During the solidification process of the alloy, the fine grains formed on the surface of the mold will fall off, and the solidified dendrite arms will also be remelted and fractured, resulting in a large number of free equiaxed dendrites in the liquid phase area, which will move under the action of natural convection and gravity [1]. The movement and falling process of a large number of group grains not only has an important influence on the formation of positive segregation at the top of ingot, A-type and V-type segregation [2], but also is the main reason for the formation of triangular cone-shaped negative segregation at the bottom of large ingot [3,4].Therefore, it is of great significance to add the calculation of grain movement process into the numerical model of ingot macrosegregation to improve the prediction accuracy.
At present, the phase field method is mostly used to simulate grain movement in the world. In 2008, Do-Quang M et al. [5] used the phase field-virtual domain method to simulate the growth and movement of single dendrite under the action of gravity. In 2012, karagadde and Bhattacharya [6] used the enthalpy method (EF) to calculate the growth of dendrites, the volume of fluid method (VOF) to calculate the movement behavior of dendrites, and the immersion boundary method (IBM) to deal with the solid-liquid interface. The simulation results show that the multi dendrite growth pattern is significantly different from that of the single dendrite. In 2013, Medvedev et al. [7] used Phase Field-Lattice Boltzmann method (LBM) coupling model to calculate the dendrite growth and

CA Model
The solute equilibrium model (ZS model) proposed by Zhu [14] is used in this paper. The growth driving force is the difference between the equilibrium crystallization concentration and the actual liquid concentration at the interface. The actual liquid concentration C L can be calculated by LBM, and the equilibrium crystallization concentration C eq L can be calculated by the following formula according to the equilibrium crystallization theory: (1) where C 0 is the initial concentration of the alloy; m is the slope of the liquidus; T L is the actual temperature of the interface; T eq L is the liquidus temperature at the initial concentration of C 0 , Γ is the Gibbs Thomson coefficient; K is the average curvature at the solid/liquid interface, f (φ, θ 0 ) is the anisotropic function of the interface energy. K can be calculated from the spatial distribution of the interface solid phase ratio.
According to the Gibbs-Thomson formula, the interface energy anisotropy function in Equation (1) can be expressed as: where δ = 15ε is the anisotropy coefficient (ε is the anisotropic strength of the interface energy), ψ(φ, θ 0 ) is the anisotropy function of the interface energy, φ is the angle between the normal direction of the solid-liquid interface and the horizontal direction, and θ 0 is the preferred growth direction. The anisotropy function ψ(φ, θ 0 ) and the growth angle φ can be calculated from Equations (4) and (5).
When C L < C eq L , the solid fraction increment ∆ f s in a time step is calculated by the following formula: In order to partially eliminate anisotropy, this paper uses an improved eight-neighbor capture method proposed by Zhu [15].

LBM Model
According to Boussinesq's approximation, the effect of latent heat and solutes on the density during the solidification process can be expressed by the following formula: where ρ 0 , T 0 , and C 0 represent the initial density, temperature, and concentration of the liquid phase, respectively, and T and C are the temperature and concentration of the liquid phase at the current Crystals 2020, 10, 70 4 of 15 moment. β T and β C are the volume expansion coefficients of temperature and concentration changes, respectively. The resultant force of the fluid particles is: The distribution function of the flow field can be expressed as: where F i is the component force of the particle under the external force field in the i direction, and its magnitude is expressed as: where ω i is the weight coefficient in each direction, which represents the probability of particles moving in different directions, which can be expressed as: The equilibrium distribution function and relaxation time of the flow field in Equation (9) are respectively expressed as: τ f = 3ν/(c 2 ∆t) + 0.5 (13) where ν is the dynamic viscosity of the fluid. The macroscopic density ρ and velocity u are obtained by adding the distribution function.
The distribution functions of temperature field and solute field are similar to those of flow field: where h i (x, t) and g i (x, t) are the distribution functions of temperature field and solute field at the position x at time t, respectively, and h i (x, t) and g i (x, t) represent the temperature field and the equilibrium distribution function of the solute field, respectively defined as: Crystals 2020, 10, 70 The relaxation time τ α in the temperature field and the relaxation time τ D in the solute field can be obtained by using the corresponding diffusion coefficients: where α is a temperature diffusion coefficient, and D is a concentration diffusion coefficient. Macro temperature and concentration are: The source term H i of the temperature field and the source term G i of the concentration field can be expressed as: In the formula, ∆T and ∆H respectively represent the latent heat released by the solidification of the alloy and the excluded solutes.

Ladd Method to Calculate the Solid-liquid Interface Interaction Force
In Figure 1, x b is the particle boundary point, x l is the liquid phase lattice point of the boundary node along the c −i direction, x s is the solid phase lattice point of the boundary node along the c i direction, and u b is the particle velocity. The calculation formulas of x s and x b are: where ∆t is the time step, V b and W b are the translational and rotational speeds, respectively, and x c is the center of mass of the solid particles. The solid-liquid distribution functions at time t are f i (x l , t) and f −i (x s , t). After ∆t/2 time, the two particles move to the boundary and collide. The distribution function at this time is: i g x t represent the temperature field and the equilibrium distribution function of the solute field, respectively defined as: e u u u x (18) ) 2 3 ) ( 2 The relaxation time α τ in the temperature field and the relaxation time D τ in the solute field can be obtained by using the corresponding diffusion coefficients: where α is a temperature diffusion coefficient, and D is a concentration diffusion coefficient.
Macro temperature and concentration are: The source term i H of the temperature field and the source term i G of the concentration field can be expressed as: In the formula, T Δ and H Δ respectively represent the latent heat released by the solidification of the alloy and the excluded solutes. In Figure 1, b

Ladd Method to Calculate the Solid-liquid Interface Interaction Force
x is the particle boundary point, l x is the liquid phase lattice point of the boundary node along the i c − direction, s x is the solid phase lattice point of the boundary node along the i c direction, and b u is the particle velocity. The calculation formulas of s x and b x are: After the ∆t/2 time, the fluid particles bounce to the corresponding lattice points respectively. At this time, the distribution functions of the liquid phase and solid phase lattice points are: Crystals 2020, 10, 70 6 of 15 Force exerted by fluid particles on solid particles: The total force F on the solid particles is: The force moment on the solid particles is: According to Newton's second law, the translation speed and rotation speed of the grains can be calculated respectively as: where M S and I S are mass and moment of inertia, respectively, and G is the combined force of gravity and buoyancy.

Processing of Solute Fields at Moving Boundaries
In order to accurately calculate the solute field during movement, a solute extrapolation method for calculating the solute field at the moving boundary is proposed.
As shown in Figure 2, the white and black grid points are the liquid and solid grid points, the solid ellipse is the position of the dendrite at the previous moment, the dotted ellipse is the current position. The gray grid points are the covered liquid nodes because of dendrite movement. The concentration C L of the covered liquid grid points will be distributed in a certain proportion, which is: As shown in Figure 2, the white and black grid points are the liquid and solid grid points, the solid ellipse is the position of the dendrite at the previous moment, the dotted ellipse is the current position. The gray grid points are the covered liquid nodes because of dendrite movement. The concentration L C of the covered liquid grid points will be distributed in a certain proportion, which is: respectively represent the concentration assigned to three lattice points along the direction of dendrite movement. The distance between the lattice points is the grain  L ,∆C 2 L , and ∆C 3 L respectively represent the concentration assigned to three lattice points along the direction of dendrite movement. The distance between the lattice points is the grain movement distance. A, B, and C are the distribution coefficients and satisfy the following relations: In this paper, the undetermined coefficient method is used, and the distribution coefficients of different proportions are used to calculate. Then the results are compared with the experimental results of Liu's single dendrite drop. The empirical values of A, B, and C are 0.7, 0.2, and 0.1, respectively. The assigned concentration will be added to the three lattice points along the direction of dendrite movement.

Verification
First, the accuracy of the Ladd method for processing moving boundaries is evaluated through the settlement process of circular particles in an infinitely long pipe. Then, the accuracy of the solute distribution model established in this paper is verified by calculating the solute field conservation of the single dendrite in the moving state. The LBM model, the dendrite growth of CA model, and the rationality of the single dendrite movement model have been confirmed in relevant researches [16,17], so this article will not be described here.

Settling of a Circular Particle in an Infinitely Long Tube
The circular particle accelerates to settle down in the infinite tube due to gravity. As the speed increases, the resistance of the ball increases. When the three forces of resistance, gravity, and buoyancy are balanced, the particle falls at a uniform speed. Glowinski pointed out in the literature [18] that, for a circular particle settling in an infinitely long pipe, when its physical parameters are determined, the settling speed is also determined at a steady state.
When the liquid is at a low Reynolds number, the resistance Ff to circular particles moving in an infinitely long pipeline is directly proportional to the settlement speed [19] and can be expressed as: where v is the falling velocity of circular particles, η is the dynamic viscosity of the fluid, K is the correction factor, and its value reflects the influence of the resistance of the pipe wall facing the particles. K is expressed as: where W i = W/D, W is the width of a long square tube and D is the diameter of the circular particle. When particles are in a state of three forces equilibrium, the resistance can be obtained as: ρ S and ρ L are solid density and liquid density, and g is gravity acceleration. According to Equations (37)-(39), the final falling speed of the circular particle can be calculated as: The size of the circular pipe in this paper is 4 cm × 8 cm. The side length of each grid is 0.01 cm, with a total of 400 × 800 grids. The particle diameter is taken as 0.48 cm, the fluid density is ρ L = 1.0 g/cm 3 and the solid density is ρ S = 1.02 g/cm 3 , and the fluid viscosity is η = 0.33 g/(cm·s). The relaxation time τ is taken as 0.8, and the gravity acceleration is g = 980.0 cm/s 2 . At the initial moment, the circular particle are placed at points (2 cm, 6 cm), and they are at rest before the calculation starts.
The simulation results are shown in Figure 3. The particle dropped due to gravity, and two symmetrical vortices formed on both sides of the particle, which is consistent with what Do-Quang [5] Crystals 2020, 10, 70 8 of 15 described in his literature. According to the simulation value of particle settlement speed, it can be seen that the simulation results in this paper agree well with the theoretical analytical solution. Therefore, the method used in this paper can be used to calculate the moving boundary problem.
Equations (37)-(39), the final falling speed of the circular particle can be calculated as: The size of the circular pipe in this paper is 4 cm × 8 cm. The side length of each grid is 0.01 cm, with a total of 400 × 800 grids. The particle diameter is taken as 0.48 cm, the fluid density is The simulation results are shown in Figure 3. The particle dropped due to gravity, and two symmetrical vortices formed on both sides of the particle, which is consistent with what Do-Quang [5] described in his literature. According to the simulation value of particle settlement speed, it can be seen that the simulation results in this paper agree well with the theoretical analytical solution. Therefore, the method used in this paper can be used to calculate the moving boundary problem. Figure 4 shows the morphology of dendrite growth under different solute field treatment methods. Figure 4a is the solute distribution method and Figure 4b is the solute extrapolation method.

Calculation of the Solute Field
In the first method, the lower end of the dendrite grew faster, because the dendrite covered the original high-concentration solute domain during the settlement of the dendrite. Through the uniform distribution of solutes, the solute gradient at the upstream end of the dendrite is reduced, so it has a faster growth rate. In contrast, dendrite growth using the solute extrapolation method was more uniform. This article compares the two methods to calculate the average regional  Figure 4 shows the morphology of dendrite growth under different solute field treatment methods. Figure 4a is the solute distribution method and Figure 4b is the solute extrapolation method.

Calculation of the Solute Field
In the first method, the lower end of the dendrite grew faster, because the dendrite covered the original high-concentration solute domain during the settlement of the dendrite. Through the uniform distribution of solutes, the solute gradient at the upstream end of the dendrite is reduced, so it has a faster growth rate. In contrast, dendrite growth using the solute extrapolation method was more uniform. This article compares the two methods to calculate the average regional concentration at the time of solidification time of 0.3 s, 0.5 s, 0.7 s. It is found that the solute field of Method 1 changes greatly, and the average solute concentrations at the corresponding moments are 4.682%, 4.673%, 4.658% (standard is 4.7%), solute field conservation is always 4.7% in Method 2. It can be seen that the solute field treatment in the second method is more accurate.

Multiple Dendrites Rotation
At the initial time, 5 crystal nuclei with preferential growth angles of 0 and a solid phase rate of 0.2 were placed in the middle of the simulation region. When the number of solidification grids is 3000, the dendrites stop growing. At this time, a uniform rotating flow is applied in the simulation region. It can be seen from the Figure 5 that the dendrite can maintain its original shape even after being rotated for multiple turns.

Multiple Dendrites Rotation
At the initial time, 5 crystal nuclei with preferential growth angles of 0 and a solid phase rate of 0.2 were placed in the middle of the simulation region. When the number of solidification grids is 3000, the dendrites stop growing. At this time, a uniform rotating flow is applied in the simulation region. It can be seen from the Figure 5 that the dendrite can maintain its original shape even after being rotated for multiple turns.

Multiple Dendrites Rotation
At the initial time, 5 crystal nuclei with preferential growth angles of 0 and a solid phase rate of 0.2 were placed in the middle of the simulation region. When the number of solidification grids is 3000, the dendrites stop growing. At this time, a uniform rotating flow is applied in the simulation region. It can be seen from the Figure 5  This article gives the formula for calculating the solid phase rate error, as follows: where 0 S f is the initial solid phase rate of the dendrite, Sin f is the solid phase rate after dendrite rotation. According to Equation (41), it is calculated that the change of the solid phase rate of the equiaxed crystal is maintained within 0.2% after 10 revolutions. In the calculations of this paper, the amplitude of dendrite rotation is small (rad <1, rad represents the radian of dendrite rotation), so it can be considered that the numerical model for calculating dendrite rotation established in this paper has no great influence on the dendrite morphology.

Discussion
This article is based on the following assumptions. The melt is an incompressible Newtonian fluid. The solid and liquid phases have the same thermal conductivity, and the heat transfer This article gives the formula for calculating the solid phase rate error, as follows: where f S0 is the initial solid phase rate of the dendrite, f Sin is the solid phase rate after dendrite rotation. According to Equation (41), it is calculated that the change of the solid phase rate of the equiaxed crystal is maintained within 0.2% after 10 revolutions. In the calculations of this paper, the amplitude of dendrite rotation is small (rad <1, rad represents the radian of dendrite rotation), so it can be considered that the numerical model for calculating dendrite rotation established in this paper has no great influence on the dendrite morphology.

Discussion
This article is based on the following assumptions. The melt is an incompressible Newtonian fluid. The solid and liquid phases have the same thermal conductivity, and the heat transfer boundary conditions are adiabatic boundary conditions. Solute diffusion only occurs in the liquid phase, there is no solute diffusion in the solid phase, and the boundary condition for mass transfer is the non-diffusion boundary condition [13].

Single Dendrite Movement
As shown in Figure 6, the dendrite arm at the front end of the dendrite grows faster in the moving state, and a secondary dendrite arm is generated. The upper dendrite arm has a lower speed and almost stops growing. The reason is that the dendrite movement compresses the lower solute boundary layer, the concentration gradient becomes larger, and the temperature gradient becomes lower; the upper solute boundary layer is stretched, the concentration gradient becomes smaller, and the temperature gradient increases. Therefore, the growth of the dendrite arms at the upstream side will be further promoted. As the falling speed increases, this asymmetric growth phenomenon becomes more and more obvious. In actual solidification, when the dendrite grows asymmetrically, the dendrite will rotate. From Figure 6c 2 , it can be seen that compared with pure falling dendrites (Figure 6b 2 ), the former dendrite arms have prominent tips, asymmetric growth is more obvious, and the latter dendrites are more uniform, which is consistent with the literature [6].
The solute distribution is also closely related to the flow of the solution around the dendrites. It can be seen from Figure 7 that the flow direction of the solution is different in the three movement states. The solution near the stationary dendrite flows upward, and the solution far away from the dendrite flows downward. For a moving dendrite, due to the high viscosity of the solution, a downward pulling force will be generated on the surrounding solution during the drop of the dendrite. This causes the solution near the dendrite to flow downwards and the solution away from the dendrites to flow upwards. In addition, for moving dendrites, two vortices will be generated behind it.  As shown in Figure 6, the dendrite arm at the front end of the dendrite grows faster in the moving state, and a secondary dendrite arm is generated. The upper dendrite arm has a lower speed and almost stops growing. The reason is that the dendrite movement compresses the lower solute boundary layer, the concentration gradient becomes larger, and the temperature gradient becomes lower; the upper solute boundary layer is stretched, the concentration gradient becomes smaller, and the temperature gradient increases. Therefore, the growth of the dendrite arms at the upstream side will be further promoted. As the falling speed increases, this asymmetric growth phenomenon becomes more and more obvious. In actual solidification, when the dendrite grows asymmetrically, the dendrite will rotate. From Figure 6c2, it can be seen that compared with pure falling dendrites (Figure 6b2), the former dendrite arms have prominent tips, asymmetric growth is more obvious, and the latter dendrites are more uniform, which is consistent with the literature [6].
The solute distribution is also closely related to the flow of the solution around the dendrites. It can be seen from Figure 7 that the flow direction of the solution is different in the three movement states. The solution near the stationary dendrite flows upward, and the solution far away from the dendrite flows downward. For a moving dendrite, due to the high viscosity of the solution, a downward pulling force will be generated on the surrounding solution during the drop of the dendrite. This causes the solution near the dendrite to flow downwards and the solution away from the dendrites to flow upwards. In addition, for moving dendrites, two vortices will be generated behind it.    It can be seen from Figure 8 that the solution flows faster around the moving dendrites. When the solidification time is 0.5 s, the flow velocity of the solution around the dendrite in the stationary, pure falling, and rotating falling states is 0.008 mm/s, −0.01 mm/s, −0.009 mm/s. This is because the rotation of the dendrite makes the contact surface larger between the lower end and the solution, resulting in an increase in the resistance of dendrites, which reduces the falling speed of the rotating dendrite to a certain extent.

Multi-Dendrite Movement
As the dendrites grow, the solute boundary layers around the dendrites come into contact and fuse with each other to form a high-concentration solute domain. The growth of dendrite arms in this domain is inhibited, and the growth of dendrite arms away from this domain is promoted. The asymmetry of the dendrite is aggravated, and the rotation speed of the dendrite is increased.

Multi-Dendrite Movement
As the dendrites grow, the solute boundary layers around the dendrites come into contact and fuse with each other to form a high-concentration solute domain. The growth of dendrite arms in this domain is inhibited, and the growth of dendrite arms away from this domain is promoted. The asymmetry of the dendrite is aggravated, and the rotation speed of the dendrite is increased.
The rotation of dendrites is different from that of single dendrites during solidification. See Table 2 for the change of preferential dendrite growth angle. The initial preferred growth angle of No. 4 dendrite is 0.1 rad. The solutes between the dendrites are close to each other and form a high concentration solute domain, therefore the growth of the right dendrite arm of No. 4 dendrite is inhibited, and the growth of the left dendrite arm is promoted. This causes the mass on the left side of the dendrite to be greater than the mass on the right side, so the dendrite is subjected to a counterclockwise torque. When the solidification time is 0.4 s, the dendrite rotates 0.3 rad in a counterclockwise direction. For the No. 1 dendrite, the lower part of the simulated domain is less affected by the solute field. It starts to rotate at 0.2 s and falls to the bottom at 0.5 s. At this time, the preferred growth direction of the dendrite is −0.5 rad. For the No. 9 dendrite, the dendrite arm at the lower end of the dendrite extends to the high-concentration domain, and the growth is inhibited. Therefore, the dendrite is more symmetrical on the left and right, and the torque is less. During its movement, the dendrite only rotated 0.2 rad. The No. 6 dendrite is located in the center of the simulation area. The uncertainty of its growth behavior leads to the complexity of its movement behavior. It can be seen that the dendrite rotates 0.3 rad clockwise when the solidification time is 0.4 s. At 0.5 s, the dendrite rotates counterclockwise by 0.1 rad, and the dendrite rotates left and right.
Multi-dendritic effects are also manifested in interactions between fluids. As shown in Figure 9a, as the dendrites grow, the fluid vortices begin to merge, forming a strong convection between the dendrites. This strong convection hinders the vertical drop of the dendrite and has a lateral force on the dendrite. The dendrite starts to move laterally, and the grains on both sides have a tendency of centrifugal movement. And the movement of the grains in the central domain shows a trend of swinging left and right. counterclockwise torque. When the solidification time is 0.4 s, the dendrite rotates 0.3 rad in a counterclockwise direction. For the No. 1 dendrite, the lower part of the simulated domain is less affected by the solute field. It starts to rotate at 0.2 s and falls to the bottom at 0.5 s. At this time, the preferred growth direction of the dendrite is −0.5 rad. For the No. 9 dendrite, the dendrite arm at the lower end of the dendrite extends to the high-concentration domain, and the growth is inhibited. Therefore, the dendrite is more symmetrical on the left and right, and the torque is less. During its movement, the dendrite only rotated 0.2 rad. The No. 6 dendrite is located in the center of the simulation area. The uncertainty of its growth behavior leads to the complexity of its movement behavior. It can be seen that the dendrite rotates 0.3 rad clockwise when the solidification time is 0.4 s. At 0.5 s, the dendrite rotates counterclockwise by 0.1 rad, and the dendrite rotates left and right. Multi-dendritic effects are also manifested in interactions between fluids. As shown in Figure 9a, as the dendrites grow, the fluid vortices begin to merge, forming a strong convection between the dendrites. This strong convection hinders the vertical drop of the dendrite and has a lateral force on the dendrite. The dendrite starts to move laterally, and the grains on both sides have a tendency of centrifugal movement. And the movement of the grains in the central domain shows a trend of swinging left and right.  Figure 10 shows the falling speed of the four dendrites. Each dendrite undergoes a process of acceleration and deceleration, and their absolute speeds are lower than those of the single dendrites. It is shown that during the growth of multi-dendritic, due to the interaction between the dendrites, a part of the gravity is offset and the dendrite's moving speed is reduced to a certain extent.  Figure 10 shows the falling speed of the four dendrites. Each dendrite undergoes a process of acceleration and deceleration, and their absolute speeds are lower than those of the single dendrites. It is shown that during the growth of multi-dendritic, due to the interaction between the dendrites, a part of the gravity is offset and the dendrite's moving speed is reduced to a certain extent.
Compared with the growth of multiple dendrites in the moving state, the growth mode of the dendrites in the stationary state is much simpler. The growth of the dendrites in the middle of the simulated domain is suppressed, and the growth of dendrites near the boundary is promoted. (As shown in Figure 11) Compared with the growth of multiple dendrites in the moving state, the growth mode of the dendrites in the stationary state is much simpler. The growth of the dendrites in the middle of the simulated domain is suppressed, and the growth of dendrites near the boundary is promoted. (As shown in Figure 11) acceleration and deceleration, and their absolute speeds are lower than those of the single dendrites. It is shown that during the growth of multi-dendritic, due to the interaction between the dendrites, a part of the gravity is offset and the dendrite's moving speed is reduced to a certain extent.
Compared with the growth of multiple dendrites in the moving state, the growth mode of the dendrites in the stationary state is much simpler. The growth of the dendrites in the middle of the simulated domain is suppressed, and the growth of dendrites near the boundary is promoted. (As shown in Figure 11) From the above analysis, the movement behavior of multi-dendrites is related to the melt convection and solute overlap between the dendrites. In response to this phenomenon, scholars proposed to apply an external force field (such as an electromagnetic field) to the dendrite to offset gravity, thereby changing the movement state of the dendrite, and then controlling the solute segregation of the casting [20,21]. From the above analysis, the movement behavior of multi-dendrites is related to the melt convection and solute overlap between the dendrites. In response to this phenomenon, scholars proposed to apply an external force field (such as an electromagnetic field) to the dendrite to offset gravity, thereby changing the movement state of the dendrite, and then controlling the solute segregation of the casting [20,21].

Conclusions
In this work, a CA-LBM-Ladd coupling model for calculating multi-dendritic motion was established, and an extrapolated distribution method for calculating the solute distribution around the dendrite under the state of dendritic motion in the melt was proposed, which realizes the calculation of the solute field in the real state of dendrite movement. The CA-LBM-Ladd coupling model was verified, and then the motion of single and multiple dendrites was simulated using this model. The simulation results show that: 1) the superposition of the flow fields between the multiple dendrites causes the movement state to change; 2) the superposition of solute field results in the change of concentration gradient, changes the growth mode of dendrite, and then changes the movement state. This is quite different from the growth mode of single dendrite.

Abbreviations
A detailed list of symbols abbreviations used in this paper is listed and a brief explanation is given.

Symbols
Unit Meaning f s mass% Solid fraction C L mass% Actual solute concentration at interface C eq L mass% Liquid-phase equilibrium crystallization concentration at the interface ∆ f s mass% Solid fraction increment k -Equilibrium partition coefficient of solute C 0 mass% Initial concentration of the alloy T L K Actual temperature of the interface T