## 3 Validation of Numerical Models
### Numerical Simulation of Three-Point Bending of Ice Beams  
Based on [PERSON] et al. 's experiment [48], a numerical simulation of the three-point bending of ice beams is carried out. Before we begin modeling, we must make two assumptions. On the one hand, the fixed support and indenter of the ice beam three-point bending experiment are metal materials with a substantially higher stiffness than ice. Therefore, they are modeled as rigid bodies. This treatment ensures that the fixed support and indenter are not deformed during the simulation, and this assumption has a negligible impact on simulation accuracy. On the other hand, there will be initial defects within the actual ice body, such as microcracks and air bubbles, and the presence of these defects is ignored during the modeling process with the assumption that the ice body has no initial defects. However, studies such as wing de-icing in the later stages of this paper are concerned with the macroscopic fracture characteristics of the ice body; the initial defects of the ice body have been avoided to the most significant degree possible when preparing the ice beam samples. Hence, this assumption has little effect on the accuracy of the simulation results. Figure 3 displays the sample size, numerical model, and experimental setup. Two fixed rigid bodies support the ice beam, while the other rigid body is loaded downward at a speed of 0.08 mm/min at the center of the upper surface of the beam. [PERSON] et al. [49] conducted three friction studies between stainless steel pin planes and ice of varying roughness. The results showed that the friction coefficients of the ice and the stainless steel were 0.02-0.11. The friction coefficient between the ice beam and the stiff body, where both are in universal contact, is set to 0.045 in this research. C3D8R is used to mesh the ice beam in the modeling process, and zero-thickness cohesive elements are introduced in the same manner described in Section 2.2.1. Table 1 displays the material properties of the zero-thickness cohesive element and the C3D8R. Figure 3: Numerical modeling of three-point bending of ice beams. This work models ice beams using three different element sizes: 4 mm, 4.5 mm, and 5 mm. It is advised that the number of elements in the cohesive zone is between 2 and 10, depending on the sample size [52; 53]. [PERSON] et al. [54] claim that to characterize the material's fracture process more precisely, there must be at least four cohesive elements inside a cohesive zone length. In order to satisfy the accuracy requirements, the cohesive zone in the ice beam modeling has to contain at least six cohesive elements. Figure 4 compares the experimental and numerical results of the ice beam's three-point bending load-deflection relationship, demonstrating that the numerical load-deflection curve for the ice beam with an element size of 4 mm is mainly consistent with the experimental results. The numerical peak load is 141 N, with a compression deflection of 0.17 mm. The relative inaccuracy is 5.53%. As the element size decreases, the slope of the load-deflection curve and the peak load decrease. Similar to the results of [PERSON] et al. [55], in Figure 5, the stress cloud is displayed. As seen in Figure 5,b, compressive stresses were created in the specimen's upper layer and tensile stresses in its lower layer during the three-point bending loading, according to Figure 5, which corresponds with the experiment's center fracture and further supports the model's validity, when loading increases, the cohesive element at the ice beam's middle bottom first experiences damage before being destroyed.