Abstract
Taking the surface profile of C45 steel grinding as the research object, a two-dimensional finite element contact model of rough surface and smooth rigid surface was established by using ANSYS software, and the stress distribution characteristics and rules of contact surface were analyzed under 0–10 MPa normal load. On this basis, the force-magnetic coupling model was established by using ANSYS APDL language. Furthermore, the influence of stress under on leakage magnetic field of contact surface under different load conditions was studied. The results indicated that according to the zero-crossing point of the normal component of the leakage magnetic field or the extreme point of the tangential component, the number of stress concentrations on the contact surface and the stress level of the corresponding area can be effectively determined. This innovative approach offers valuable insights for future studies on surface contact stress distribution between components.
Introduction
In mechanical assemblies, the surface contact between parts was not smooth surface-to-surface contact, but consists of mutual contacts between many micro-convex bodies [1–3]. The real contact area will be much smaller than the nominal Contact area. When the equipment was running, the contact surfaces between parts were prone to vibration, reducing the stiffness of the equipment. According to research, the stiffness of the joint surface in a machine tool accounted for approximately 60% to 80% of the total stiffness of the machine tool, and the deformation caused by the joint surface accounts for approximately 85% to 90% of the total static deformation of the machine tool [4]. When two surfaces squeezed each other due to the action of load, the actual contact area will be subject to a large contact load, and phenomena such as contact surface wear and local crushing were prone to occur [5]. Therefore, conveniently and accurately obtaining reliable surface contact stress distribution is a key theoretical and technical issue in the design and development of general mechanical products, and it is also a research hotspot in domestic and foreign academic circles.
Under certain assumptions, Fu et al. revealed the nonlinear relationship between the normal contact stiffness of the joint surface and the parameters of the joint surface based on the contact fractal theory [6]. Wang et al. proposed that the wavelet decomposition method can be used to identify the contour characteristic length scale parameters in the fractal parameters of rough surfaces, and verified the effectiveness by comparing with other methods [7]. Zhao et al. proposed a new elastic-plastic micro contact model of rough surfaces based on the theory of contact mechanics to study the elasticity, elastic-plasticity, plastic contact deformation, contact area and contact load of asperities on rough surfaces [8]. Wang et al. established a fractal contact stiffness model considering the interaction of asperities based on the proposed concept of fractal smoothness, and derived the relationship between the deformation of a given asperity and the contact stiffness [9]. In recent years, magnetic non-destructive testing technology based on the magnetization theory of ferromagnetic materials has developed rapidly [10]. This technology can not only quickly detect the stress state of materials, but also evaluate the damage state of materials [11–13]. Roskosz et al. applied the magnetic memory method to evaluate the stress distribution of specimens containing holes under tension [14]. The research results showed that there is a clear correlation between the residual magnetic field and the stress state. Lopez et al. studied the impact of the stamping process on residual stress and magnetic domain structure [15]. The results showed that the residual compressive stress caused by the stamping process will change the shape and width of the magnetic domain structure near the edge of the tested sample. Yang et al. used the high magnetic permeability characteristics of ferromagnetic materials to establish an evaluation system for measuring defect characteristics through magnetic flux leakage signals [16].
Numerical simulation was an effective method to analyze the quantitative relationship between leakage field distribution and damage or other physical properties of ferromagnetic materials, and was an important part of theoretical research on ferromagnetic materials [17]. Sorabh et al. used ANSYS finite element software to simulate the leakage flux distribution on the surface of the pipe after the residual stress and large concave deformation of the ferromagnetic pipe under the action of spherical pressure load [18]. Ren et al. simulated the force-magnetic coupling process of 18CrNi4A steel and obtained the change rules of the magnetic memory signal under different loads [19]. Li et al. used numerical simulation software to analyze the relationship between defects in pipeline steel and leakage flux density, which provided a basis for quantitative analysis of pipeline steel defects [20]. Shen et al. used Magnet software to conduct finite element static simulation of a pipeline magnetic flux leakage detector and studied the effects of defect size and sensor lift-off value on the magnetic flux leakage signal [21]. Numerical simulation methods can more intuitively describe the magnetic field and stress distribution on the contact surface through finite element models, and can more accurately present the simulation results of magnetization [22–24].
Our research goal was to use the force-magnetic effect to characterize the stress distribution characteristics of two rough surfaces under a certain load. We assumed that the stress distribution level in the contact area can be judged based on the magnetic field distribution pattern of the asperities on the rough contact surface. Using ANSYS finite element software, a plane model of the rough surface was first constructed, and the contact properties between the rough surface and the rigid plane under normal load (≤10 MPa) were studied. On this basis, ANSYS APDL language was used to establish a force-magnetic coupling model to further simulate the impact of different load conditions on the leakage magnetic field on the contact surface.

Rough surface contour curve.

Finite element contact model. (a) Overall mesh of the contact model. (b) Local mesh of the contact model.
Construction of finite element contact model
Firstly, the surface profile of grinding sample of C45 steel was obtained by TR300 surface profile measuring instrument. The diameter of the probe of the surface profile measuring instrument is 2 μm, and the measurement is scanned along the surface perpendicular to the grinding texture. Figure 1 showed the surface profile of the grinding sample of 45 steel. The sampling interval was 1 μm, the sampling length was 5 mm, and a total of 5000 discrete sampling point height values were obtained. The surface roughness (R a ) of the sample was 1.20 μm, and contour root mean square deviation (R q ) was 0.93 μm. The measured data points were processed using the least squares method and connected with straight lines to form a two-dimensional rough surface profile. Further, on the basis of two-dimensional surface profile data, APDL language was used to program the profile data into ANSYS, and the finite element plane strain contact model of two-dimensional rough surface was established, as shown in Figure. 2. In order to reduce the number of meshes as much as possible, simplify the calculation and shorten the calculation time, a two-dimensional finite element contact model of rough surface and smooth rigid surface was established. the model gradually densified the mesh in layers along the vertical direction, making the mesh density transition appropriately. In order to ensure the calculation accuracy and the convergence of the contact analysis, a uniform quadrilateral grid was used near the contact surface to evenly divide the surface contour in the horizontal direction, ensuring the uniformity of the grid in the horizontal direction, as shown in Figure 2(b). The height and width of the model were 5 mm, which was consistent with the length of the sampled rough surface contour data. The number of units in the entire model was 38647, the number of nodes was 96024, and the unit type was PLANE183. The material properties of the model were shown in Table 1, using a multilinear elastic-plastic strain hardening material model obtained from uniaxial tensile experiments. Add vertical constraints to the bottom of the model, fully constrain the lower left corner, and add horizontal constraints to the rigid surface. Finally, the loading condition is to apply a normal pressure of 10 MPa on the rigid surface of the sample, and load according to 100 loads with equal step length.
Material properties of the unit.
Material properties of the unit.
In the geomagnetic field environment, when the ferromagnetic material had no external load, the magnetic induction line was concentrated in the ferromagnet, almost no magnetic induction line comes out, and the magnetic flux was mostly parallel to the surface of the ferromagnet. When there was external load, the local stress concentration caused the change of the internal magnetic property parameter (permeability) of the ferromagnetic material, thus causing the change of the surface leakage magnetic field [14]. Based on this, according to the functional relation (1) of the permeability and stress of ferromagnetic materials under uniform magnetic field established in literature [25], the variation law of magnetic leakage signal with stress was analyzed.
Therefore, based on the finite element contact model, we obtained the stress distribution results of the model by applying different loads in the static field. Then, the permeability of ferromagnetic materials under different stress levels was obtained according to the force-magnetic relation (2). Finally, the magnetic properties of the ferromagnetic material element were assigned by APDL program in the static magnetic field, so that the surface leakage magnetic field distribution of the ferromagnetic material under the loading condition of the model can be obtained.

Magnetostatic analysis model. (a) Geometric model. (b) Grid model. (c) Enlarged grid model. (d) Enlarged mesh model of the contact area.
In order to analyze the distribution law of magnetic flux leakage signal with stress in the stress concentration area, the dimensional model adopted the same mesh unit division in the statics and magnetostatic analysis. Figure 3 showed the geometric model and mesh division of two-dimensional magnetic field. In magnetostatic analysis, an air model should also be added to the outside of the model to collect the magnetic leakage signal in the air on the outer surface where the stress was concentrated. The air mesh in contact with the contact surface was denser, and the air mesh became sparser farther away from the contact surface. The model used the PLANE53 unit, which had 8 nodes and each node has 4 degrees of freedom. The unit node output included magnetic field intensity component, magnetic field intensity vector, magnetic flux density component, and magnetic flux density vector. In order to obtain the magnetic leakage signal size of each node of the unit, we chose the magnetic field intensity vector (Az) as the output result. The boundary conditions of the model were as follows: the upper air was the geomagnetic field intensity Az = 39.8 A/m, and the lower air was the lower limit Az = 0 A/m [26].
Contact surface analysis of the model
Figure 4 showed the stress cloud diagram of 100 load steps. From the stress cloud diagram, it can be observed that the stress distribution of the contact layer, especially the surface stress close to the contact point, was very complex. The stress started from the contact point and gradually decreased inward. The accelerated splitting near the contact layer itself caused freezing caused by the contact, and also included the extrusion and puffing caused by other contact points.

Stress cloud. (a) Overall stress cloud after 100 loading steps. (b) Stress cloud diagram of area A.

Relationship between displacement and normal load.
Figure 5 showed the relationship curve between normal deformation and normal load, in which the abscissa was the normal load value and the ordinate is the displacement in the normal direction of the midpoint of the rigid surface. It can be seen from the figure that the deformation displacement gradually increases with the increase of normal load, and the relationship between the two was generally linear. L.H. He established the finite element contact model between rough surface and ideal rigid surface, and compared the finite element calculation results with the experimental results to verify the rationality of the model [4]. By comparing with the experiment in reference [4], the test result curves of the two were consistent, indicating that the model calculation results were reliable.
Figure 6 showed the relationship curve between contact area and normal load, where the abscissa was the size of the normal load and the ordinate was the contact area (contact line length). It can be seen from the figure that when the load was less than 1 MPa, the contact area increased generally linearly with the increase of the normal load. When the load exceeded 1 MPa, the increase in the contact area became smaller and smaller, approaching saturation.

Relationship between contact area and normal load.
Figure 7 was the result curve of magnetostatic analysis, in which figures (a)–(c) were respectively the tangential and normal leakage magnetic field intensity curves of air 1 mm above the rough surface under various loading conditions. In Figure 7(a), the abscissa was the coordinate value in the actual length direction, and the ordinate was the tangential leakage magnetic field intensity of the air layer above the corresponding position. In Figure 7(b), the abscissa was the coordinate value in the actual length direction, and the ordinate was the normal leakage magnetic field intensity of the magnetic field in the air layer above the corresponding position. Figure 7(c) synthesized the leakage magnetic field under normal loads from 1.0 MPa to 10 MPa. Figure 7(d) combined the leakage magnetic field under various normal loads from 1.0 MPa to 10 MPa, and the differences in the leakage magnetic field under various load conditions can be clearly compared. From Figure 7, we found that the stress value in the contact area was greater than the stress value in the uncontact area, resulting in the magnetic permeability of the contact area being smaller than the uncontact area, which caused magnetic leakage in the contact area, thereby increasing the normal component of the magnetic field intensity in the leakage magnetic field. All cross zero, and the tangential components all have the maximum value of the curve shape. As the load increases, each curve in the figure had more tangential magnetic field components showing maximum values, and more normal components showing zero crossing points. This is because more contact areas appear as the load increases.

Distribution of leakage magnetic field above the contact surface under different loads. (a) Tangential leakage magnetic field strength under 0.1–1.0 MPa load. (b) Normal leakage magnetic field strength under 0.1–1.0 MPa load. (c) Tangential leakage magnetic field strength under 1.0–10 MPa load. (d) Normal leakage magnetic field strength under 1.0–10 MPa load.
As shown in Figure 8, the stress cloud diagram of the contact area under 10 MPa normal load, it can be observed that there were 4 obvious contact areas, which were related to the number of extreme points of the tangential leakage magnetic field intensity and the normal leakage magnetic field intensity under 10 MPa. The number of zero points was the same.

Stress cloud in contact area under normal load of 10 MPa.
As shown in Figure 9, under the normal load of 10 MPa, the leakage magnetic field intensity of air layered with different heights from the contact surface. It can be seen from the figure that as the height increases, the values of the tangential and normal leakage magnetic field strengths of the air layer above the contact surface gradually decrease, which indicated that the lifting height has an impact on the leakage magnetic field signal.

Leakage field strength. (a) Tangential leakage field strength. (b) Normal leakage field strength.
By using the force-magnetic coupling model established by ANSYS, we can judge the number of stress concentrations in the contact area, the real stress distribution and the stress in the contact area based on the zero-crossing points of the normal component and the extreme points of the tangential component of the leakage magnetic field size, which was consistent with our hypothesis. At the same time, the research showed that the amplitude of tangential and normal leakage magnetic field intensity of the air layer above the contact surface decreases gradually with the increase of height. Finally, we effectively proved the applicability of the magnetic permeability effect of stress and the force-magnetic coupling model, as well as the feasibility of studying the surface contact characteristics of parts.
