Multi-Physics Finite Element Model of Relay Contact Resistance and Temperature Rise Considering Multi-Scale and 3D Fractal Surface

Contact resistance of the main circuit is an important parameter for measuring the performance and reliability of relay products. The contact resistance value interacts with the stress and temperature of the contact surface during the relay’s operation, and it will change with the ambient temperature and load current. Therefore, it is difficult to calculate the contact resistance when the relay is in the working environment. However, the finite element model used to predict contact resistance in previous studies does not consider the effect of current and temperature on the contact surface stress at the microscopic level. In this paper, based on the fractal theory, we established three-dimensional electrical contact finite element models in micrometer-scale and nano-scale to solve the contradiction between the computational efficiency and accuracy of the finite element model containing fractal surfaces. In the microscopic electrical contact model, the multi-physics coupling process of electric-temperature-stress at the conductive spot is analyzed, to improve the accuracy of the calculation of the local contact resistance at the conductive spot. Through the data transfer between the macroscopic finite element temperature rise model of the relay and the microscopic electrical contact model, the temperature rise and contact resistance of the relay under different load currents are calculated more accurately. The feasibility of the simulation method is proved by comparison with experimental measurement results.


I. INTRODUCTION
The relay is a commonly used switching device, and its main circuit contact resistance is an important parameter. High-power relays are increasingly used in electric vehicles and other fields. As the size of the relay decreases and the power increases, the temperature rise caused by the contact heating becomes more and more obvious. Especially when used in extreme temperature environments such as aerospace, the high temperature of the relay in the working state will affect the coil resistance and ferromagnetic material properties. This in turn affects the electromotive force of the electromagnetic coil, significantly affects the working characteristics of the device, and even causes it to fail or reduce its service life. The contact resistance value of the relay directly affects the accuracy of the heating calculation. The contact The associate editor coordinating the review of this manuscript and approving it for publication was Giambattista Gruosso . resistance value is related to the contact's own material and processing technology, and interacts with various factors such as stress and temperature, making the calculation of the contact resistance of the relay in work more difficult. Contact resistance is affected by ambient temperature and load and needs to be predicted in advance. However, the internal temperature rise of the relay and the contact resistance of the relay in different working environments are difficult to measure in real time. When measuring in the laboratory, it also faces the problems of a long cycle and high cost of making prototypes and redesigning. Therefore, it is necessary to calculate the contact resistance and simulate the temperature field of the relay in the working state.
There have been many studies on contact resistance. Earlier Holm proposed a basic model of contact resistance. Based on the statistical description of rough surfaces, Greenwood established a G-W model based on the Holm model, and pointed out that there were both elastic deformation VOLUME 8, 2020 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ and plastic deformation during the contact process of micro convex spots, and the concept of plastic index was introduced [1], [2]. Mandelbrot proposed the concept of fractal, and then proposed W-M fractal function based on Weierstrass function [3]. The micro-profile of the contact surface conforms to the fractal feature. Majumdar and Bhushan proposed the M-B model based on the W-M fractal function to simulate the rough contact surface, and simplified the model to a model where the rough surface contacts a smooth rigid plane [4]. However, in the calculation of the deformation amount δ, the compression process of the micro convex spot was not considered. Wang and Konmvopoulos proposed a modified model of the M-B model, which gave the relationship between the contact area and the contact cross-sectional area of micro-contacts [5], [6], which is widely used later. However, these models use 2D contours to represent rough surfaces and are only applicable to isotropic surfaces. Yan and Komvopoulos established a 3D model of fractal contact surfaces for contact analysis of irregular surfaces [7]. In Morag's model, as the contact pressure increases, the deformation process of the micro convex spot was transformed from elastic deformation to plastic deformation, which is divided into three stages: elastic contact, elastoplastic contact and plastic contact [8]. This type of elastoplastic contact problem usually equates the shape of the contact spot to a circle. Kogut analyzes the effect of thin insulation film and current on contact resistance between rough surfaces [9]. L. Pei et al. Established a three-dimensional finite element model of the fractal surface to analyze the relationship between the actual contact area and the load [10]. The fractal method began to be used more in the calculation of electrical contact resistance [11], [12] and thermal contact resistance [13]- [15]. Wilson used stacked elastic-plastic spheres to establish a multi-scale model of contact area and load relationship [16]. Capelli et al. used genetic algorithm to optimize the calculation of fractal parameters of rough surfaces when using fractal method to calculate constriction resistance [17]. Wanbin Ren et al. adopted the analytical method and the finite element method to simulate the electrical contact resistance of coated spherical contacts, but did not consider the micro level of the contact surface [18]. Dong et al. established a surface layer contact model considering rough surfaces, but this layered structure was slightly different from the contact surface of the relay contacts [19], [20]. Wahab et al. established a finite element model of the contact position, analyzed the stress and wear at the micro-scale contact point, and used the multi-scale method to improve the calculation efficiency [21]- [25]. Hong Liu and J. W. McBride measured the micro-topography profile of the contact surface and established a finite element model of a spherical micro-contact with a rough surface, and analyzed the change of the contact area when it comes in contact with a plane with the loading and unloading of the contact pressure [26]- [30]. However, the contour data obtained through the measurement are affected by the accuracy of the machine's sampling. Using the data directly for modeling may affect the fineness of the model. The fractal method is more suitable for modeling the micro-surface of relay contacts with random micro-surface contours and large apparent contact area. In addition, the above-mentioned research work did not consider the process that the contact resistance value interacts with factors such as temperature and stress when the load current was switched on. Based on fractal theory and contact surface topography measurement data, a 3D finite element model of the irregular topography of the contact surface at the microscale was established. Fractal parameters are obtained by repeatedly measuring the contours of different positions on the contact surface, which can make the modeling of irregular surfaces as accurate as possible. A microscopic electrical contact finite element model including irregular contact surface and rigid plane was constructed. The multi-physics coupling process of electric-temperature-stress in the micro-scale finite element model of electrical contact was considered when the contacts were energized. According to the self-similar property of the surface, as the sampling scale decreases, the contact model will appear as numerous spot contacts with zero contact area. When the sampling length is large, with the increase of the fineness of the finite element model, too many mesh units will make the model difficult to solve. In this study, a method of constructing a multi-scale model was introduced to solve the problem that it was difficult for the micro-scale finite element model of electric contact to take into account both calculation accuracy and calculation efficiency. Then, the contact resistance of the relay contacts under different load currents is calculated through the data exchange between the micro electrical contact model and the macro relay temperature rise finite element model. The contact resistance of the relay was calculated as an interpolation function of the temperature rise of the contact surface according to the contact surface morphology and the current level being switched on. This function was used to calculate the overall temperature rise of the relay during steady state operation. Finally, the comparison and analysis between simulation data and experimental data was carried out.

II. MICRO-SCALE FRACTAL SURFACE MODEL
The macro-scale contact model of the relay contact is a plane-to-plane contact, but the micro-scale contact form is shown in Figure 1. Current passes through multiple conductive spots, and the actual contact resistance is the parallel connection of the contact resistances of multiple conductive spots. However, the conductive spots are not completely in contact, and are still a parallel connection of the contact resistance of multiple conductive spots on a more microscopic level. It is difficult to construct a fine finite element model containing a large number of tiny contact spots, and the model is difficult to solve, especially when multi-physics coupling is considered. Therefore, two kinds of micro-scale surface finite element models were established in this paper, namely micrometer-scale and nano-scale, both of which are constructed using fractal functions. Copper moving contact of the relay used in the experiment and the contact surface measured by 3D digital microscope and the contour curve on a vertical cross section on the surface are shown in Figure 2. In this relay, due to different processing techniques, the surface of the moving contact is rougher than the static contact. Therefore, the surface of the static contact is considered as a rigid plane, and the surface of the moving contact is regarded as a rough surface. The middle part of Figure 2 shows a 3D microscopic image of the surface of moving contact at a magnification of 1000 times.

A. FRACTAL FUNCTION
Fluctuations of micro-scale contact surface contour are highly random, and the surface topography data of a single observation may not be a good representation of the characteristics of the entire contact surface. The W-M function can be used to describe the fractal surface contour with continuity and self-affinity [4].
where Z (x) is the height of the rough surface profile, G is the fractal roughness, D is the fractal dimension (1< D <2), and γ n is the spatial frequency (γ > 1), which determines the frequency spectrum of the surface roughness. When the height distribution of the rough surface meets the Gaussian distribution, take γ = 1.5. In addition, n m is the ordinal number of the lowest cutoff frequency, with γ n m = 1/L, where L is the sampling length.

B. CALCULATION METHOD OF FRACTAL PARAMETER
Fractal parameters can be calculated by structural function method, power spectrum method and other methods. The structure function method was used in this study. The measured surface contour curve is regarded as a time series Z (x). A time series with fractal characteristics can make the structure function of the sampled data fit: where τ represents an arbitrary selected value of the data x interval, and S (τ ) represents the arithmetic mean of the variance. Derived from the self-affinity of W-M function and (1), where the coefficient C is Draw a double logarithmic plot of S (τ ) and τ . If it is a straight line, get the slope k s of the straight line and the intercept B of the straight line on the log S (τ ) axis.
According to the formulas (5) and (6), D and G are obtained. The three-dimensional fractal dimension can be approximately converted from the two-dimensional fractal dimension calculated by D s = D + 1.

C. 3D FRACTAL SURFACE
This paper uses the three-dimensional W-M function proposed by Yan and Komvopoulos [7] to build a 3D model of the rough surface.
where D is the three-dimensional surface fractal dimension (2 < D < 3), M is the number of fold overlaps of the fractal surface, and φ m,n are Gaussian random numbers. The other parameters are the same as in (1). After repeatedly measuring the contact surface contours at different sampling positions and calculating by the structure function method, the fractal parameter is D = 2.3, G = 1.28 × 10 −9 m. when L = 4 × 10 −4 m, the 3D rough surface constructed by the fractal function is shown in Figure 3. The roughness of the fractal surface generated by the W-M function is close to the measured surface. There are irregular tiny pits on the surface of the measured contact due to the processing technology, most of which are located in the wave valleys of the surface, which have little effect on the actual contact area.

III. MULTI-PHYSICS AND MULTI-SCALE MODEL OF RELAY CONTACT RESISTANCE AND TEMPERATURE RISE
Generally, contact resistance is considered to be composed of two parts: constriction resistance and film resistance. The relay used in this experiment has the following characteristics: 1) The contact is in a closed arc extinguishing chamber and does not contact with oxygen, so it is not easy to form an oxide film on the contact surface.
2) The relay has been switched many times in advance, which can damage the oxide film on the surface of the contacts to a certain extent.
3) The contact passes a large current, which has a damaging effect on the oxide film.
Based on the above considerations, the film resistance has little effect on the contact resistance of the relay contact, so only the constriction resistance is considered in this article.
The multi-scale model in this paper includes macro-scale, micrometer-scale and nano-scale. The macro-scale is used for the macroscopic temperature rise model of the relay, and the micrometer scale and nano-scale are used for the modeling of irregular micro-scale surfaces.

A. MACRO-SCALE FINITE ELEMENT MODEL OF RELAY TEMPERATURE RISE
The macro-scale model of relay temperature rise is a millimeter scale model, and its finite element model and boundary conditions are shown in Figure 4. The relay is lifted by a thick wire and does not contact other objects. This can simplify the boundary conditions of the model and facilitate the comparison between experimental data and simulation data. The ambient temperature is controlled at 293K. Cooling conditions include natural convection and radiation. The formula is as follows.
where is the heat flux, h is the convection heat transfer coefficient, A is the heat dissipation area, ε is the material emissivity, and σ is the Stefan-Boltzmann constant with a value of 5.67 × 10 −8 .
The heat source includes the main circuit Joule heat and the coil Joule heat. The main circuit heating is the Joule heat when the current flows through the moving contact, the contact resistance of contacts, the static contact and the wire terminal. The setting method of the contact resistance on the contact surface of the contact is to set the boundary condition formula: where h c is the constriction conductivity, which indicates the conductivity of the contact surface, and its unit is S/m 2 . J is the current density. V 1 and V 2 are the potentials of the contact surfaces on both sides. The relationship between the contact resistances R c and h c is: where A c is the apparent contact area between the moving contact and the static contact. The model is easy to solve but the structure is complex, so the element type used is free tetrahedral. Each element has four nodes. When the other conditions of the model are the same and only the element sizes of the model are different, the temperature simulation results at point C1 are as follows. When the number of model elements is 106218, 560886, and 3156935, the temperature at the C1 point is 95.996K, 96.370K, and 96.476K, respectively. The element size in the macro model has little effect on the simulation results. The minimum element size used in the model is 0.452mm, and the total number of elements is 560886. The main purpose of establishing the macroscopic temperature rise model of the relay is to calculate the temperature field of the relay and the temperature T u of the contact surface of contacts when the relay works until the temperature rise reaches a steady state.

B. MULTI-PHYSICS FINITE ELEMENT MODEL AT MICRO SCALE
Micro-scale finite element models include micrometer-scale models and nano-scale models. The two models are basically the same except for the size. To take into account both calculation accuracy and solution efficiency, it is necessary to select a suitable model size, that is, the sampling length L for fractal surface modeling. In order to make the calculation results accurate, the 3D finite element model of the microscopic surface of the contact requires fine meshing, which will generate a large number of meshing units. The contact resistance of the entire contact surface is a parallel connection of the contact resistances of multiple single microscale models. When the sampling length is large, a single microscale model is large in size and can contain more contact spots, but the model is difficult to solve. When the sampling length is small, the single model is small in size and contains a small number of contact spots. As the number of models in parallel increases, the total number of contact spots increases. This makes some points at the troughs in the surface contour contact the rigid plane in the model, which is not in line with the actual contact situation. Referring to the surface morphology observed under the 3D digital microscope, L u = 400µm was selected as the sampling length of the micrometer-scale contact model in this study. L n = 400nm was selected as the sampling length of the nano-scale contact model.
Taking the micrometer-scale model as an example, its finite element model and boundary conditions are shown in Figure 5. The model is equivalent to a smooth rigid plane in contact with a 3D fractal rough surface [4]. The elastoplastic deformation of the rough plane is considered. The Young's modulus E on the rough plane is given by where E is the composite Young's modulus [2], E 1 and E 2 are the Young's modulus of the materials on both sides, and υ 1 and υ 2 are the Poisson's ratios of the materials on both sides. The micrometer-scale model represents a small part of the contact surface of the macro model, as shown in Figure 7(a). The temperature at the conductive spot of the micromodel is higher than the rest part. Because the height of the model is small enough, the temperature of the upper and lower surfaces of the micrometer-scale model is set to the temperature T u of the contact surface in the macro model. For the nanoscale model, the temperature of the upper and lower surfaces of the model is set to the temperature T n of the contact spot  of the micrometer-scale model. Adjacent to the micro-scale model is the same model, in which temperature and current are hardly transmitted laterally. Therefore, the surrounding surfaces are set to electrical insulation and thermal insulation.
The total contact pressure F m of the relay contacts, and the apparent contact area A c of the dynamic and static contacts are all constant values. The current I m is the total current value passed by the contacts. The side length of the micrometerscale model, that is, the sampling length is L u = 400µm. It is assumed that the contact surfaces of the moving contact and the static contact are parallel when in contact. The current through the microscale model I u = I m L 2 u /A c . An electrical contact pair is provided between the rigid plane and the fractal surface, so that current can pass through the part where the two surfaces are in contact with each other. Set the electric potential of the lower surface of the model to 0 and the current I u flow into the upper surface to calculate the voltage between the upper and lower surfaces. In order to make the model easy to converge, the following setting method is adopted. Surface to surface contact pairs are set on the two contact surfaces in the model. Fix the upper and lower surfaces of the model so that the two contact surfaces squeeze each other, resulting in elastoplastic deformation and vertical contact force F c . The length d p in Figure 5 represents the distance between the upper and lower surfaces of the model, which represents the distance between the upper half of the model The parameter transfer between different physics in this model is as follows. The model heats up after passing the current, resulting in a temperature rise T . The effect of temperature on contact resistance is divided into two aspects, namely the effect on the resistivity of the material and the effect on the elastoplastic mechanical parameters of the material. The effect of temperature on copper resistivity ρ is referred to the formula: where ρ 0 is the material resistivity when the temperature rise is 0, α is the temperature coefficient of resistance, and the α value of copper is 0.00393/K. As the temperature increases, the interatomic distance increases, and the Young's modulus of the metal material gradually decreases. When the temperature does not reach 0.52T m , the Young's modulus of the metal decreases approximately linearly with increasing temperature. When the temperature rise is T , its Young's modulus Where E 0 is the Young's modulus of the material at a temperature rise of 0, and σ is the thermal expansion coefficient of the material. When the temperature is 293K, the Young's modulus of copper E 0 = 1.1 × 10 11 Pa, and the thermal expansion coefficient is σ = 17.5×10 −6 /K. T m is the melting point of metallic materials, and the melting point of copper is 1357K.
The material of the moving contact is hard oxygen-free copper. The material parameters at the temperature of 293K are shown in Table 1. The temperature also affects the yield stress σ s of the material, and its curve with temperature is shown in Figure 6.
The change in Young's modulus E T and the yield stress σ s of the material after the temperature rise will affect the elastoplastic deformation of the rough surface, which will increase the area of the contact spots. In addition, the change in metal conductivity will affect the constriction resistance of the entire contact model. Because the current input to the model is constant, the change in constriction resistance is reflected in the potential difference between the upper and lower surfaces of the model. According to the heating power Q = I 2 R of the contact spot, the change in resistance will cause the change of the model temperature. The temperature at the contact spots of the model will be higher than the temperature at the upper and lower surfaces. The nanoscale model is the same as the micrometer-scale model, and both involve the coupling of electro-thermal-structural physical fields. The micrometer-scale model mainly provides calculation results such as the constriction coefficient h c , the temperature T n at the contact spot, and the contact pressure F n at the contact spot. The calculation results provided by the nano-scale model are mainly the contraction conductivity h cn at the contact spot of the micrometer-scale model.
The model also uses free tetrahedral meshing, with the smallest element size on the fractal surface. The element size has a significant effect on the calculation results. Because the irregular fractal surface needs a fine mesh to construct. When other simulation conditions are the same, the comparison of contact force simulation results of models with different element sizes is shown in Table 2. As the element size decreases, the accuracy of the simulation results increases. When the element size of the contact position on the fractal surface reaches 7.2µm, the simulation results already have a large deviation. Using larger element sizes will make the simulation results more inaccurate, and at the same time non-convergence begins to appear. When the minimum element size is 0.6µm and 1.6µm, the simulation results are relatively close. Smaller element size means more solution time. In order to take into account the computational efficiency, the smallest element size we use is 1.6µm.

C. DATA TRANSFER BETWEEN MULTI-SCALE MODELS
The parameter transfer process between macro-scale model, micrometer-scale model and nano-scale model is shown in Figure 7(a). The locations of heat sources in different scale models are also marked in blue in the figure.
The data transfer process between multi-scale models is shown in Figure 7(b). First, set the initial values of material properties and boundary conditions for each scale model. The ambient temperature in the simulation is 293K. The constriction conductivity h cn calculated from the nano-scale model is transferred to the micrometer-scale model. h cn is set on the contact surface of the contact spot of the micrometer-scale model to calculate the potential difference U u between the upper and lower surfaces of the model. Then the constriction conductivity h c for the apparent contact surface of moving contacts and the static contacts is derived. Then transfer h c to the macro-scale model for calculating the Joule heat of contacts. The heat generated by the coil and the load circuit in the macro -scale model increases the temperature of the model, and the temperature T u of the contact surface in the model is transferred to the micrometer-scale model as its temperature boundary condition. The temperature T n and pressure F n at the contact spot are calculated from the micrometer-scale model and transferred to the nano-scale model as its temperature and force boundary conditions. The new constriction conductivity h cn calculated from the nano-model is transferred to the micrometer-scale model. This completes a cycle of data exchange. Iterate continuously until the model converges to get the calculation results.
However, the multi-scale model is complicated and contains too many elements, which is difficult to converge and has extremely low calculation efficiency. Therefore, the interpolation function method is adopted in this study to improve the calculation efficiency. First, the parameter h cn calculated from the nano-scale model is used as an interpolation function of the temperature T n and the contact pressure Then h c can be substituted into the macro-scale model as an interpolation function of T u and I m . When the current I m is input, the value of h c at steady state and the temperature field of the relay can be obtained through iterative calculation. According to (11), the contact resistance value of the contact surface can be obtained from the h c value.
In this method, key parameters such as h c are transformed into an interpolation function table by a method of solving a micro-scale model in a parametric scanning manner. By calling the interpolation function, the frequent data transfer process between models is replaced, which greatly reduces the calculation cost. By setting h c = f (T u , I m ) as the boundary condition at the heat source of the macro temperature rise model, the temperature field and contact resistance when the relay works to a steady state can be quickly calculated directly in the macro-scale model.

IV. COMPARISON OF SIMULATION AND EXPERIMENTAL RESULTS
A certain type of relay used in this experiment is shown in Figure 8. When the main circuit is connected to a constant current, the voltages at points A and B are measured. The resistance between the two points can be obtained by R = U /I . The experimental ambient temperature was 293K. In the experiment, 4 relays of the same type were measured to exclude errors caused by accidental factors. When the temperature rise is 0, the resistance values between the A and B points of the four relays are 0.136m , 0.136m , 0.138m , and 0.142m . There are two static contacts and one moving contact in each relay, so there are two contact surfaces in the circuit between the moving contact and the static contact. VOLUME 8, 2020 Excluding the contact resistance, the conductor resistance of the copper circuit consisting of moving contacts and static contacts is 0.022m . The contact resistances of the two contact surfaces in the load circuit are in series, so the contact resistances of each contact surface of the four relays are 0.057m , 0.057m , 0.058m and 0.06m . According to the method of this paper, the contact resistance of a single contact surface is 0.0493m when the temperature is 297K and the current is not turned on. The calculation error is 13.5%-17.8%.
Because the contact is in a sealed arc-extinguishing chamber, its temperature cannot be measured, so the temperature of the contact surface of the contact is calculated by simulation. The temperature outside the relay is measured with a thermocouple. The four sampling positions of the thermocouple probe are marked as C1-C4 in Figure 8. The point C1 is outside the static contact and is close to the voltage sampling point A. Point C2 and point C3 are on the side of the relay case. Point C4 is below the center of the bottom of the relay.
The cloud diagram of the temperature rise simulation results of the macro-scale finite element model of the relay when the load current is 200A is shown in Figure 9. The sampling position of the temperature T u of the contact surface is also marked in the figure. When the relay is connected to different load currents, the steady-state temperature of the four thermocouple sampling points is shown in Figure 10. The simulation results of the macro-scale model of relay temperature rise are similar to the measured results, which proves that the model can predict the temperature of the contact point well. In addition, it can be known from Figure 9 that due to the good thermal conductivity of copper, the temperature outside the contact is close to the temperature at the internal contact surface.
As shown in Figure 10, as the current increases, the error between the simulation and the actual measurement increases. There are three reasons for the error. 1. The probes of the thermocouples are all attached to the outside of the relay. Although thermal grease has been used, the measurement results are slightly lower. 2. In this study, the contact thermal resistance between the internal components of the relay was not considered, resulting in a slightly smaller total thermal resistance of the heat transfer path from the contact to the relay case. As a result, the simulated temperature of the case is slightly higher, and the simulated temperature of the contacts is slightly lower. 3. The boundary condition parameters are calculated by the formula, and there will be some errors compared with the actual situation. The contact is the heat source of the relay, and its heating power is P = I 2 R. When the current increases, the power of contact heating increases, and the errors caused by points 2 and 3 also increase with the increase of current. Considering contact thermal resistance and improving the accuracy of boundary condition parameter calculation in future research can improve the accuracy of temperature model calculations.
The simulation cloud diagram of the micrometer-scale finite element model when the contacts are connected with a load current of 100A to the steady state is shown in Figure 11. The figure includes the results of potential, temperature, and von mises stress of the micrometer-scale model.
The simulated values of contact surface temperature and the comparison between the measured results and the simulation results of the contact resistance when four relays are connected to different load currents to reach a steady state are shown in Figure 12. The point where the current is 0 in the figure (a) represents the temperature of the contact surface when the load current is not turned on and only the coil is energized to reach a steady state. The simulation method in this paper can be used to calculate the contact resistance of relay contacts that contact in a plane-to-plane form. It can be seen from the figure that as the load current increases, the steady state temperature of the contact surface increases, and the contact resistance shows a downward trend. In this study, it is assumed that the contact surfaces of the moving contact and the static contact in the model are parallel. In fact, although the contact surface of the contact has been adjusted for parallelism during the assembly process, it is difficult to achieve a completely parallel state. This makes the actual contact area slightly smaller, which is one of the reasons that the simulation value of the contact resistance is smaller than the actual value.

V. CONCLUSION
This paper presents a multi-scale and multi-physics calculation method of contact resistance by using a fractal function to establish three-dimensional micro-scale finite element models of irregular contact surfaces and combining with a macro-scale model of relay temperature rise. By changing the boundary conditions such as the ambient temperature of the model and the current value of the main circuit load, this method can be used to calculate the steady-state temperature rise and contact resistance of the contact under different temperature environments when different load currents are switched on. By comparison with experimental measurements, it can be seen that the calculation error of the contact resistance calculation method proposed in this paper is within an acceptable range and has certain application value. With the increase of current, the temperature of the relay contacts increases, and the contact resistance value decreases. By adopting a method in which the parameters that need to be transferred between models, such as the constriction conductivity h c , are made into a mutual interpolation function, the convergence rate and calculation efficiency can be greatly improved. The effect that the contact surfaces of the moving contact and the static contact are not completely parallel is not considered in the model, and the calculated contact resistance may be slightly smaller because of this effect. Therefore, the influence of assembly factors on the calculation of contact resistance can be added to future research. For the case where the contact is connected with a larger load current and a larger deformation occurs at the contact spot, the mesh-free method can be considered in the micro-scale model in future research.