A numerical solution to the effects of surface roughness on water–coal contact angle

Coal dust is a great threat to coal mine workers' health and safety in coal mine production. Wet dust removal is one of the effective dust removal methods. As a solid, coal has different rough surfaces, which have a certain effect on the wetting effect of coal. In this paper, three coal samples with different surface wettability are used as the research objects. Phase-field interface tracking method is used to simulate the wetting of droplets on rough surfaces. From the simulation results, it can be concluded that the influence of the rough interface on the contact angle of the droplets is in accordance with the change rule described in the Wenzel model. As the roughness increases, the contact angle of the hydrophilic lignite surface gradually decreases. As the roughness increases, the contact angle of hydrophobic coking coal gradually increases. The change trend of the contact on the surface of weakly hydrophilic anthracite coal is the same as that of lignite. Due to the local and global differences, the contact angles obtained from the numerical model are slightly different from the values calculated from the Wenzel model.


Scientific Reports
| (2021) 11:459 | https://doi.org/10.1038/s41598-020-80729-9 www.nature.com/scientificreports/ The study also compared the influence of gas-liquid ratio on the contact angle, and found that the gas-liquid ratio has the largest contact angle within a certain range. The effect of surface micro/nano structure on wetting process was studied through the simulation of lattice Boltzmann method 10 . It is found that when the droplet size exceeds the size of the microstructure, the contact angle can be calculated by the classic Cassie equation. When the sizes of the two are almost the same, the contact angle cannot be simply calculated by the Cassie equation. Under different mining conditions or stress states, coal presents a surface with different roughness, for example, the surface of drilling and excavation, due to the different coal quality and operation mode will show different morphology. It will affect the contact area of water on the coal surface, and thus affect the wetting effect of water on coal. Besides chemical composition of coal surface, coal surface physical property like roughness or morphology can also play an important role in wetting process. However, there is a slight lack of research on the influence of coal surface roughness on wetting, and there are fewer studies on numerical simulation. As a conventional research method, numerical simulation has particularly mature simulation methods and simulation software. When studying the roughness of the coal surface, numerical simulation can ensure the uniformity of the roughness value, and the influence of roughness on the wetting of the coal surface can be studied qualitatively and quantitatively. Therefore, the purpose of this paper is to study the effect of roughness on the wetting of the coal surface via a numerical solution.

Coal sample selection and natural fracture roughness analysis
Coal sample selection. Three presentative coal samples were selected, lignite from Hami, coking coal from Anyang, and anthracite from Zhaogu, which have hydrophilic, weakly hydrophilic and hydrophobic properties. Some important parameters of the coal sample, such as the proximate analysis and the hardness coefficient f, were tested according to the methods specified in the national standards GB/T 212-2008 and GBT 8208-1987. The results are shown in Table 1.
Roughness of natural fracture surface. The surface roughness is measured using the optical contact Angle 3D morphology coupling instrument produced by Biolin Scientific. The roughness measurement uses the 3D Topography Module, which is composed of a sample platform, a grating projector, a digital camera, and an information processor. The 3D Topography Module is based on Phase Measurement Profilometry (PMP). This technology uses a projector to project a number of sinusoidal grating projections on the surface of the measured object, and the sampling camera collects the grating projection images. Because the collected information is modulated by the measured object, it is necessary to obtain the wrapping information of the height and the position through the demodulation algorithm, and use the unwrapping algorithm to unwrap to obtain the unwrapping information of the position. Finally, through the conversion relationship between the position and the height, the height data of the object surface is obtained 11 . The process of the experiment is shown in Fig. 1. www.nature.com/scientificreports/ The morphology and roughness of the natural fracture surfaces of the three coal samples were measured on relatively flat fracture surface to ensure convenient measurement. For each coal sample, three surfaces were selected to be tested within the area range of 1.4 mm × 1.0 mm (Fixed measuring area of equipment), and each surface was scanned five times, and then use the mean square roughness value as the evaluation index to analyze the difference in surface morphology and roughness of different coal samples. Figure 2 shows the optical images of the three coal samples. Figures 3 and 4 presents the reconstructed 2D topography and 3D topography of the fractured surfaces of the three coal samples respectively. Due to the reflection of coal sample, some data is missing, and the missing place will be displayed as blank.
According to the above morphology measurement, the roughness of the natural fractures of the three coal samples are Hami lignite R q = 6.54 μm, Anyang coking coal R q = 5.14 μm, and Zhao anthracite R q = 3.52 μm (where R q is the mean square roughness). From the perspective of coal hardness coefficient, different coal samples have natural fracture surfaces with different roughness.  www.nature.com/scientificreports/ Numerical model Theoretical model. The concept of contact angle was proposed to quantitatively characterize the wettability of solid surfaces. Generally the surface with a contact angle of less than 90° is defined as a hydrophilic surface, and the surface with a contact angle of greater than 90° is deemed a hydrophobic surface 12 .
The contact angle formulation on the ideal homogeneous and smooth surface was derived by the Young's equation 13 , which relates the intrinsic contact angle to the surface tension between the three phases.

Contact angle θ:
Young's equation can be used to predict the ideal contact angle of absolute smoothness, but it is limited to characterizing the contact angle of a droplet on an ideal solid surface. The smaller the contact angle, the better the wettability, and the larger the contact angle, the worse the wettability.
The wettability of the real solid surface was studied and introduced the surface roughness coefficient, namely the 'roughness factor (λ)' 14 . Wenzel also corrected the Young's equation to establish a model of the liquid completely wetted in the local groove area, the equation is written as: or where θ r is the apparent contact angle, λ is the ratio of the actual surface area to the projected area and is always greater than 1. Therefore, in combination with Eq. (3), it is known that the absolute value of the contact angle cosine function formed by the interaction between a real solid surface and a droplet is larger than that of an ideal smooth solid surface. When the solid surface is hydrophilic (contact angle < 90°), the apparent contact angle will decrease as the solid surface roughness increases and the hydrophilicity is better. If the solid surface is hydrophobic (contact angle > 90°), the apparent contact angle will increase with the solid surface roughness. However, it was also found that the hydrophilic surface becomes hydrophobic under certain conditions and it was difficult for Wenzel to explain this phenomenon 15 .
In the study of surface wetting, it was found that during the wetting process of rough solid surfaces, the liquid did not completely enter the surface gaps, and there was still a certain amount of air under the droplets, and the droplets were suspended between the air and the solid composite contact surface 16 . They analyzed from the perspective of thermodynamics and proposed the following equation: where θ r represents the apparent contact angle of the composite contact surface, θ 1 and θ 2 represent the intrinsic contact angle of the two media, f 1 and f 2 represents the proportion of the two media at the composite interface, When the composite contact surface is composed of air and solid, since the contact angle of liquid and air is 180°, Eq. (5) can be simplified to: (3) cosθ r = cosθ. Because wetting is the interaction in the gas-liquid-solid interface, the key to the simulation is to track the three-phase (gas-liquid-solid) interface. There are currently several methods for phase interface tracking, such as Volume of Fluid (VOF) 17 , Lattice Boltzmann Method (LBM) 18 , Particle-in-cell method (PIC) 19 , Level Set method 20 , Phase Field method 21 . The phase field method has practical physical meaning, which is based on the kinetic process of thermodynamics. In order to characterize the changes of moving objects in time and space, ϕ(r, t) variable is introduced, r represents space, and t represents time. ϕ= −1 indicates the liquid phase zone, ϕ=1 indicates the solid phase zone, and a value between -1 and 1 indicates the solid-liquid two-phase zone. This method also uses the Cahn-Hilliard equation 22 to control the phase field change, which is more suitable for the calculation of microfluidic motion and has high calculation accuracy. Therefore, this paper chooses the phase field method to numerically simulate the wetting of water droplets on the coal surfaces.
Physical model. The simulation content is the process of wetting and spreading a droplet on a rough solid surface in an air atmosphere. The established model is a two-dimensional model, the calculation area is 18 mm × 6 mm, the drop is in the center, and the volume is 10 μL (close to the size of the experimental drop). The geometric model is shown in Fig. 5. The roughness surface is obtained by programming, and the random rough surface curve data points are obtained. This method can set different roughness value curve data through parameter adjustment, and then import the data into the COMSOL Multiphysic geometry to obtain the random rough surface, as shown in Fig. 6. where ρ is the density (kg/m 3 ), µ is the dynamic viscosity coefficient (Ns/m 2 ),p is the pressure (Pa),I is the identity matrix,g is the Acceleration of gravity (m/s 2 ),Fs is the surface tension of the solid-liquid interface (N).
In the phase field, the surface tension in Eq. (9) can be calculated according to the following equation: where φ is the phase field parameters, G is the chemical potential (J/m 3 ). Cahn-Hilliard equation: tracing the diffusion interface of two immiscible fluids. The diffusion interface is defined as the area where the dimensionless phase field variable φ ranges from − 1 to 1. In the simulation process, the Cahn-Hilliard equation is divided into two equations:  www.nature.com/scientificreports/ where γ is the mobility (m 3 s/kg), λ is the mixed energy density(N), ε is the interface thickness parameter(m), ψ is the phase field auxiliary variable. The relationship between mixing energy density and interface thickness and surface tension coefficient: Generally set the interface thickness parameter to ε = h c /2 , which h c is the characteristic grid size of the area passing through the interface. The mobility γ determines the time scale of Cahn-Hilliard diffusion, consider specific requirements when selecting, because this parameter needs to be small enough to ensure that the convection term is not excessively suppressed, and it needs to be large enough to keep the interface thickness constant.
In Eq. (10), φ is the phase field parameters that ranges from − 1 to 1. The volume fraction of each fluid: In the model fluid 1 is defined as water and fluid 2 is defined as air.
The density (kg/m 3 ) and viscosity (Pa/s) of the two-phase interface mixture: where w represents water, a represents air.

Contact angle control equation:
where n is the normal vector, n = ∇φ |∇φ| . This equation can determine the wettability of the interface.
Intrinsic contact angle measurement. An important parameter of the three coal samples is the intrinsic contact angle. The intrinsic contact angle is the contact angle formed by the substance itself and deionized water under ideal smooth surface conditions. Here we use the sessile drop method to measure the contact angle. A drop is placed on a prepared coal sample surface by using a micro syringe. The profile of sessile drop is photographed for contact angle measurement using a computer program. The intrinsic contact angle can be used to judge the hydrophilicity and hydrophobicity of the substance. Since it is troublesome to add chemical property parameters of solid materials in the simulation, it is simplified to replace the hydrophilic and hydrophobic properties of each coal sample with intrinsic contact. However, it is difficult to achieve an ideal smooth surface, so multiple high-precision sandpapers are used to grind the coal samples to make them close to the ideal smooth. The specific contact angles of the three coal samples are shown in Fig. 7, where the proximate intrinsic contact angle of each coal sample is 60° for Hami lignite, 95° for Anyang coking coal, and 85° for Zhaogu anthracite. www.nature.com/scientificreports/

Simulation conditions.
(1) Assumptions 1. There is no reaction between the wall and the liquid, no heat transfer, and the wall temperature is constant. 2. The fluid is an incompressible Newtonian fluid, and the fluid flow conforms to the laminar flow law. 3. The liquid has no penetration on the solid surface, only wetting and spreading on the surface.
(2) Initial conditions The ambient temperature is 298 K (keep at the same temperature as the experiment), the pressure is 101 kPa, the liquid phase is water defined as fluid 1, the gas phase is air defined as fluid 2. The viscosity and density parameters of the fluid derive from the selected material, and the surface tension coefficient is 0.072. The phase interface thickness parameter is determined according to the mesh size, which is one-half of the characteristic mesh size of the phase interface, and the mobility ratio takes the default value.
(3) Boundary conditions As shown in Fig. 8, the left side of the numerical simulation calculation area is set as the pressure inlet, the right side is set as the pressure outlet, and the pressure value is 0. The upper wall is set as a no sliding wall, and the lower wall is set as wetting wall. (4) Meshing The thickness of the mesh division will affect the accuracy and speed of the calculation. In order to improve the calculation accuracy and reduce the calculation time, the droplet boundary and the solid wetted wall are finely divided, and the other areas are divided with coarse mesh, as shown in the Fig. 9. The surface of the droplet and the wetting wall is dense, the smallest unit is 0.0009 mm, the largest unit is 0.06 mm, the maximum unit growth rate is 1.08, the curvature factor is 0.25, and the narrow area resolution is 1. Other calculation areas default to "more fine", the smallest unit is 0.024 mm, the largest unit is 0.16 mm, the maximum unit growth rate is 1.1, the curvature factor is 1, and its narrow area resolution is 1.  www.nature.com/scientificreports/

Model validation
According to the intrinsic contact angle value and the natural fracture surface roughness value obtained above, the simulation parameters in the COMSOL Multiphysic software are adjusted to simulate the wetting in three different situations. The real contact angle is measured by using the contact angle measurement module in the optical contact angle profiler. The model is validated by calibrating dynamic water droplet spreading process and by comparing contact angles till the steady state. The real and simulated droplet spreading processes on the Hami lignite surface are shown in Fig. 10. The processes are completed in a very short time and they are almost consistent during the time.
Comparison of contact angle between real experiment and simulation in the static state is given in Table 2. Roughness of the model is obtained by substituting the roughness value obtained from the topography measurement into the model. It is found that the experimental and simulated values are similar under the conditions of natural fractures. Taking into account the hysteresis in the real situation, the experimental value is smaller than the simulated value. In general, the contact angle simulated by this model is feasible and the model is validated to undertake more parametric studies.

Parametric study and discussion
Parametric study. In the numerical simulation, the structural roughness of the surface is 0.5 μm, 1 μm, 2 μm, 3 μm and 5 μm, and the intrinsic contact angles of 60°, 85° and 95° are used as the research objects to simulate hydrophilic, weakly hydrophilic and hydrophobic of the wetting process. The contact angle of each coal sample surface is simulated. Figure 11 list the numerical results under three different wettability conditions and the trend of contact angle changes, respectively. From the simulation result image, it can be seen that as the roughness increases, the contact angle of the hydrophilic lignite surface gradually decreases, from 58.7° to 44.3°. As the roughness increases, its hydrophilicity becomes stronger. The change trend of hydrophobic coking coal is opposite to the previous two. As the roughness increases, the contact angle gradually increases, and its contact angle increases from 96.5° to 101.7°, becoming more hydrophobic. The change trend of the surface contact boundary of weakly hydrophilic anthracite is the same as that of lignite. Its contact angle reduces from 84.7° to 76.3°, and the hydrophilicity is improved. From the perspective of hydrophilicity and hydrophobicity, Hami lignite has the best hydrophilicity, and the influence of roughness on its wettability contact angle is about 14°. Anyang coking coal has the worst hydrophilicity, and the influence of roughness on the contact angle is about 6°. The Zhaogu anthracite is in the middle, with a variation range of about 8°.
Discussion. As shown in Fig. 12, the simulation result is that there is no air gap between the droplet and the solid, which conforms to the Wenzel model. www.nature.com/scientificreports/ where λ is always greater than 1, which is the ratio of the actual surface area to the projected area. In the roughness parameter design, the value of is 1.05, 1.1, 1.2, 1.3, 1.5. According to the Wenzel model, the theoretical curve of the relationship between the roughness factor and the contact angle is drawn. In Fig. 13, for the convenience  www.nature.com/scientificreports/ of comparison, the curves with the intrinsic contact angles of 85° and 95° are partially zoomed in and listed as (a) and (b), respectively. The Wenzel model can theoretically explain the above simulation results. For hydrophilic surfaces, when the contact angle is less than 90°, the contact angle gradually decreases as the surface roughness increases, and the hydrophilicity becomes stronger. For hydrophobic surfaces, when the contact angle is greater than 90°, the contact angle gradually increases as the roughness increases, and the hydrophobicity becomes stronger. A weakly hydrophilic surface has a smaller range of changes than a hydrophilic surface. For example, when the intrinsic contact angle is 60°, as the roughness factor increases, the contact angle decreases faster than the intrinsic contact angle of 85° decreases with the increase in roughness. This also explains, with the increase in roughness, that the change range of lignite contact angle is larger than that of anthracite.
The theoretical value of the apparent contact angles calculated according to the Wenzel theoretical model are compared with the simulated value, as shown in Fig. 14. Due to the local and global difference, the simulated value is the contact angle measured on the local area of the droplet, and the theoretical value is the contact angle calculated based on the entire rough surface area. When the contact angle is about 90°, the simulated value is compared with the theoretical value: the absolute value of the simulated contact angle is more than the theoretical value, and the contact angle change rate in the simulation is greater than the theoretical value. For hydrophilic surface, the simulated contact angle is less than the theoretical value. The simulation change law of the influence of the rough interface on the contact angle of the droplet accords with the Wenzel model with slight difference though. www.nature.com/scientificreports/

Conclusion
It is validated that the proposed numerical model is a feasible and efficient tool to investigate the effect of roughness on the wetting process of water on a solid matter surface. The model is then employed to study the wetting process of coal surfaces with different hydrophilic and hydrophobic properties under the influence of roughness, and the following conclusions are obtained: (1) The numerical simulation of roughness has a certain feasibility for the wetting of coal surface. The process of droplet spreading is similar to the real situation, including the speed of droplet spreading and droplet morphology. The contact angle values under the same roughness are similar. While the numerical simulation is too ideal, resulting in the simulated contact angle value being larger than the actual situation. (2) The wetting condition of the coal surface affected by the roughness conforms to the Wenzel model. The contact angle of hydrophilic lignite decreases with the increase of roughness, and the hydrophilicity is better. The contact angle of hydrophobic coking coal increases with the increase of roughness, and the hydrophobicity is better. The change in contact angle of weakly hydrophilic anthracite is between the two, and the difference in change is also in the middle. (3) Due to the local and total difference, the simulated value is the contact angle measured on the droplet range, and the theoretical value is the contact angle calculated based on the entire rough surface range. The simulated value of the contact angle is slightly different from the theoretical value obtained from the Wenzel model. www.nature.com/scientificreports/