Numerical Simulation and Modeling of Mechano–Electro–Thermal Behavior of Electrical Contact Using COMSOL Multiphysics

1. Introduction

Electrical contacts are crucial circuit components present in most electric devices, such as relays, switches, breakers, connectors, breakers and integrated circuits. The theory and analysis of electrical contacts [1,2,3,4,5] is a multidisciplinary issue since it involves electrical, mechanical and thermal concepts, all aimed at understanding and improving the performance of these electrical connections used in various applications.

A contact system always consists of a pair of electric contacts and a source of contact force [6], with each contact made of conductive material, usually metal or metallic alloys. When a pair of contacts touch, electrical current passes with a certain electrical contact resistance. There are several factors [7] that influence the current flow and thus this resistance: contact area [8], shape [9], roughness [10], oxidation and coatings [11,12], mechanical loads [13], etc. The electrical contacts’ resistance should be low and stable, since its increase leads to contact failure and very high Joule heat losses [14].

The typical failure modes for electrical contacts consist of: (a) temperature increase in the contact area producing melting or quick aging, (b) elevated electrical contact resistance caused by fretting wear or oxidation, and (c) excessively small or large closing/opening forces.

Failures of electrical contacts and the causes leading to them are extensively studied in the scientific literature, based on several different analysis approaches.

First, there are theoretical analyses based on model types such as the single point (or Holm tube) contact model [15,16], the multi-point contact model [17,18], or the fractal geometry theory model [19,20,21]. A number of several other scholars used analytical models to study electrical contact behavior with authors of [22] presenting a study built on the improved elastic rough surface model, [23] introducing a new analytical model for the analysis of switching devices contact area for different voltage and current magnitudes and [24] proposing a semi-analytical model based on the single asperity concept to highlight the flow of current between two contacting surfaces.

Next, electrical contacts have been also widely studied using experimental setups and measurements, employing a combined analytical and experimental investigation [25] to describe the electrical-thermal behavior of electrical contact systems, using a low-cost tabletop indenter setup [26] to measure the in-situ electrical contact resistance, determining experimentally the adhesion force and energy at contact interfaces for conductive, isotropic, and rough surfaces [27]. High current automotive connectors were tested [28] to assess their thermal performance, experimental observations on contact resistance were summed up [29] to be used as an electrical contact guide, two distinct techniques were used during investigations [30] for different pressures, temperatures, and contact types, and trials were conducted [31] on copper/brass plug-and-sockets type of contacts.

Finally, there are the vast array of numerical methods. These started several decades ago with finite elements discretization [32,33] and solving in codes like FEAP [34,35,36] or discrete elements simulations using the LMGC90 open platform or MULTICOR software [37,38,39,40] but in the recent years the continuous advance of computer technology and the leap in performance and processing power allows to model, analyze, optimize and solve more and more complex problems of electrical contact systems. Nowadays, the numerical methods use dedicated engineering computer software, and are based on various techniques such as the boundary element method (BEM) [41,42,43], finite element method (FEM) [44,45,46] or coupled/multiphysics methods [47,48,49] to model and simulate contact behavior.

Angadi et al. [50] and other scholars reviewed and presented several computer software available to researchers for the analysis of electrical contacts and their behavior, with MSC Marc, ANSYS, ABAQUS and COMSOL being the most used ones as subsequently detailed:

(a) MSC Marc software was used in [51] to create the model of an electric contact and compute the temperature distribution of joule heating compared with experimental results. The performance of pin and socket electrical contacts was assessed in [52] using a computer model developed also in MSC Marc software and numerical simulation;

(b) ANSYS was used by Duan et al. [53] to develop a coupled mechanical–electrical–thermal model, in order to prove the relation between electric contact shrink range, contact resistance and maximum stress also considering the temperature. Authors of [54] proposed a multi-physics model developed in ANSYS to investigate the thermal behavior of an electrical contact with rough surface in two scenarios: direct and load transfer. The model developed introduced in [55] combines ANSYS with a MATLAB code to accurately predict the contact forces, thermal contact resistance, surface separation and ECR and provide an efficient tool to study electrical contact behavior from electrical, thermal and structural perspective;

(c) ABAQUS and FE-SAFE software are both used [56] to numerically investigate and assess in detail the behavior of aviation electrical contacts. Authors of [57] use an ABAQUS developed thermal-electric coupled simulation to study the contact pair insertion force and obtain the contact resistance. A analysis was conducted using a new electro-thermo-mechanical contact model developed by [58] in ABAQUS and two user-defined subroutines to show how the contact state significantly influences the local heating;

(d) Finally, the Swedish COMSOL Multiphysics software is used by Zhang et al. [59] to introduce a new fluid-solid heat transfer model coupled with contact mechanics and thermo-electric effect, in order to highlight the importance of the pressure on the failure of an electrical contact. Another COMSOL developed model is proposed in [60] to conduct research on the contact temperature characteristics in order to improve the performance of electrical contacts. Authors of [61] investigate by simulation in COMSOL of the effects of the electric current and contact force on the temperature and resistance of a spring type electrical contact.

However, regardless of the method applied or the software used all have both limitations and advantages, so at this point there is no imposed standard requiring a particular computer program or method of analysis. Taking this into consideration, given the extensive literature review pointing to the available software tools, and considering the experience and previous approaches [62,63,64,65] of the same research team in solving various engineering tasks using numeric simulation, in the present study COMSOL software [66] was chosen for the analysis and SolidWorks for the creation of the geometric model.

The purpose of the present paper is to develop a modeling approach based on Holm electrical contact theory in order to investigate the dependance between the electric potential and the temperature of the electrical contact by numerical methods using the COMSOL software and FEM. The results obtained are in line with the expectations and with past research presented in the literature review section.

The paper is structured as follows: Section 1 introduces the purpose if the paper and presents a comprehensive literature review on the studied thematic, Section 2 presents the theoretical aspects regarding electric contacts. Section 3 shows the initial geometric model development using SolidWorks software, followed by the step-by-step creation of the COMSOL model and the simulation studies run in the same software, and section 4 discusses the results obtained.

2. Theoretical Considerations about Electrical Contacts

As described in the Introduction section, an electrical contact consists of a pair of usually metallic contact elements, by the touch of which conduction in an electrical circuit is established. In practice, the touching of the two elements is achieved by pressing one element against the other by means of a force. This contact force can be produced by either springs (by compression or stretching) or bolts that mechanically join the two elements. In order to study and express contact behavior in mathematical formulae, it is necessary to know the main physical processes that occur when forming an electrical contact. One can imagine a simple experiment: two straight metal cylinders with similar cross section and very finely polished end surfaces are pressed against each other with a force F, as shown in Figure 1. If the measured resistance R_m along the length comprising the junction of elements l is compared with the resistance value calculated as R_c = ρ ^. l / A where l is the length, A is the transversal cross-section of the cylinders and ρ is the material resistivity, it is found that R_c << R_m. The explanation for this lies in the constriction of current lines of flow, and in the formation of a disturbing film on the processed surfaces.

No matter how finely the contact surfaces are polished, the actual touch is achieved only in certain areas called contact points, where the current lines of flow are bent and undergo a constriction as shown in Figure 2. Also, the surfaces of contact elements react with the atmosphere forming a film of Cu₂O on copper surfaces and Ag₂S on silver ones which also increases contact resistance.

In contacts with a large area the effective contact area is smaller as compared to the apparent touching area, because surfaces that appear smooth at real scale are actually rough and contaminated at microscopic scale [67]. The contact elements actually touch in micro-points of contact which are grouped into contact areas that are called contact spots or a-spots where, where the material is deformed. If the area of contact points A_F on which the contact force F is exerted is only a tiny fraction (0.01—0.05) of the apparent area A_a, there is a surface contact as illustrated in Figure 3 where contact is achieved through three a-spots noted 1, 2, and 3.

Each contact spot consists of micro-areas where the material is deformed either plastically, elastically or both. The dependence between the contact force F and the contact area A_F on which this force is exerted under plastic deformation conditions can be expressed using the formula of Holm [68]:

$$A_F=\frac{F}{\xi \cdot H}=n\cdot \pi \cdot a^2$$

(1)

where a is the radius of the equivalent circle for each of the n contact zones.

Equation (1) shows that the magnitude of the contact area does not depend on the contact force, but only on the hardness of the material H and coefficient ξ. This coefficient is subunitary and was introduced in the equation to mark that the hardness of the contact peaks is less than the hardness H measured macroscopically by pressing a sphere against a flat surface. If the touch between the two elements is achieved through a single contact area, a single-point contact is obtained. In practice this can be achieved by means of an element including a spherical cap and an element including a flat surface, as shown in Figure 4. In this case the dependence between the contact area π . a² and the contact force F is:

$$A_F=\pi \cdot a^2=\frac{F}{H}$$

(2)

Equation (2) does not include a dependence between the contact area and the radius of the contact element in the field of plastic deformation where technical contacts belong.

For the calculation of the constriction resistance two models can be applied [69], namely, the infinite conductivity sphere model and the flattened ellipsoid model.

2.1. The Model of the Infinite Conductivity Sphere

The model of the infinite conductivity sphere shown in Figure 5 consists of two half-spaces 1 and 2 of finite conductivity which model the contact elements, where the electrical conduction is established by means of a sphere of radius a and infinite conductivity. The current lines are radial, and the equipotential surfaces are spherical. The current density is constant over the surface of a sphere of radius r.

If $\rho$ is the resistivity of the material of the half-spaces, the elementary resistance is:

$$\mathrm{d}{R}_s=\frac{\rho \mathrm{d}r}{2\pi r^2}$$

(3)

and by integration the constriction resistance of one half-space becomes:

$$R_s=\frac{\rho }{\pi }{\displaystyle \underset{a}{\overset{\infty }{\int }}\frac{\mathrm{d}r}{r^2}=\frac{\rho }{2\pi a}}$$

(4)

The total contact constriction resistance of the two half-spaces in touch with the infinite conductivity sphere having a diameter 2a, is twice the value given by Euation (4), i.e.,:

$$R=\frac{\rho }{\pi a}$$

(5)

2.2. The Flattened Ellipsoid Model

In the case of the flattened ellipsoid model, the half-spaces 1 and 2 of finite conductivity are touching by means of a flattened ellipsoid, as shown in Figure 6. The equipotential surfaces are confocal ellipsoids, with the flattened ellipsoid as their base ellipsoid. To calculate the electrical resistance between the equipotential surface of the base ellipsoid (as a place of contact) and the surface of another equipotential confocal ellipsoid, formal analogy is used, which exists between the formulae characterizing the stationary electric field of direct current in a conducting medium and the formulae characterizing the electric field in an uncharged dielectric.

In this case, the expression of striction conductance is identical to the expression of capacitance in which permittivity ε is replaced by conductivity σ. Solving the problem in this way is possible assuming the same boundary conditions are allowed for both the conductive medium with conductivity σ and the dielectric medium with permittivity ε. Thus, both the base ellipsoid and the medium delimiting the surface of the confocal ellipsoid are assumed to have infinite conductivity, and the conductive medium between the surfaces of the ellipsoids has constant conductivity σ. Only in this case are the surfaces of the two ellipsoids equipotential surfaces in the stationary electric field, and the current lines are orthogonal to the surfaces of the ellipsoids.

If it is considered that the base ellipsoid is flattened as in Figure 6, and that the contact surface is an ellipse, then the equipotential surfaces in the contact members are semi-ellipsoids with the equation:

$$\frac{x^2}{{\alpha }^2+\mu }+\frac{y^2}{{\beta }^2+\mu }+\frac{z^2}{\mu }=1$$

(6)

where α and β are the semi-axes of the flattened base ellipsoid in the plane xy. Confocal ellipsoids are defined by the parameter μ (μ=0.1 ^. a², μ=0.5 ^. a², μ=2 ^. a², μ=5 ^. a², ...). Thus, the semi-axes are: $\sqrt{{\alpha }^2+\mu }$on axis x, $\sqrt{{\beta }^2+\mu }$on axis y, and$\sqrt{\mu }$on axis z.

The capacitance between the base ellipsoid and a confocal ellipsoid is:

$$C=\frac{8\pi \epsilon }{{\displaystyle \underset{0}{\overset{\mu }{\int }}\frac{\mathrm{d}\mu }{\sqrt{\left({\alpha }^2+\mu \right)\left({\beta }^2+\mu \right)\mu }}}}$$

(7)

and the capacitance between the base semi-ellipsoid and a confocal semi-ellipsoid is:

$$C_s=\frac{4\pi \epsilon }{{\displaystyle \underset{0}{\overset{\mu }{\int }}\frac{\mathrm{d}\mu }{\sqrt{\left({\alpha }^2+\mu \right)\left({\beta }^2+\mu \right)\mu }}}}$$

(8)

Based on the afore mentioned analogy, the constriction resistance between the base semi-ellipsoid and a random confocal semi-ellipsoid is:

$$R_s=\frac{1}{4\pi \sigma }{\displaystyle \underset{0}{\overset{\mu }{\int }}\frac{\mathrm{d}\mu }{\sqrt{\left({\alpha }^2+\mu \right)\left({\beta }^2+\mu \right)\mu }}}$$

(9)

If it is admitted that the contact ellipse is actually a circle, meaning that α = β = a then the elliptical integral becomes a simple transcendent integral, and:

$$R_s=\frac{1}{4\pi \sigma }{\displaystyle \underset{0}{\overset{\mu }{\int }}\frac{\mathrm{d}\mu }{\left(a^2+\mu \right)\sqrt{\mu }}}$$

(10)

With the variable changes$\sqrt{\mu }=z$,$\mu =z^2$, and $d\mu =2zdz$ we obtain:

$$R_s=\frac{1}{2\pi \sigma }{\displaystyle \underset{0}{\overset{\sqrt{\mu }}{\int }}\frac{\mathrm{d}z}{a^2+z^2}}=\frac{1}{2\pi \sigma a^2}{\displaystyle \underset{0}{\overset{\sqrt{\mu }}{\int }}\frac{\mathrm{d}z}{1+{\left(\frac{z}{a}\right)}^2}}$$

(11)

With a new variable change $\frac{z}{a}=u$, $dz=a\cdot du$ we obtain:

$$R_s=\frac{1}{2\pi \sigma a}{\displaystyle \underset{0}{\overset{\frac{\sqrt{\mu }}{a}}{\int }}\frac{\mathrm{d}u}{1+u^2}=\frac{\rho }{2\pi a}\arctan \left(\frac{\sqrt{\mu }}{a}\right)}$$

(12)

If instead of an ellipsoid a semi-space is considered then $\mu \to \infty$, and Equation (12) becomes:

$$R_s=\frac{\rho }{4a}$$

(13)

The total constriction resistance is twice as high as$R_s$:

$$R=\frac{\rho }{2a}$$

(14)

According to Holm [67], the current density across a flattened ellipsoid contact area can be approximated as:

$$J\left(r\right)=\frac{i}{2\pi a}\cdot \frac{1}{\sqrt{a^2-r^2}}$$

(15)

The analysis of the structure of electrical contacts which have bonded due to the current passing through them, shows that soldering and therefore melting of the material occurs at the edge of the point contact which justifies the acceptance of the flattened ellipsoid as a calculation model. In Figure 5 of the infinite conductivity sphere, the current lines of flow for 0.25 I, 0.5 I, 0.75 I, and I were plotted. It is found that at μ=0.5 ^. a² the confocal ellipsoid comes very close to a sphere, and at μ=20 ^. a² it is practically similar to a sphere. For this situation, the distance to the flattened ellipsoid is $z=\sqrt{\mu }=4.472\cdot a$, and the semi-axis after the axis x is $\sqrt{a^2+\mu }=4.582\cdot a$.

Knowledge of the dependence between the contact force and electrical resistance provides a first criterion for the dimensioning of electrical contacts. For a point contact in the case of elastic deformation, the contact area between two spheres or two cylindrical bars is calculated using Hertz’s formula:

$$a=\sqrt[3]{\frac{3}{4}\cdot F\cdot \left(\frac{1-{\sigma }_1^2}{E_1}+\frac{1-{\sigma }_2^2}{E_2}\right)\cdot {\left(\frac{1}{r_1}+\frac{1}{r_2}\right)}^{-1}}$$

(16)

where the notations (according to Figure 7a,c) are: a – radius of circular contact surface; F – the contact force; σ₁, σ₂ –Poisson’s ratio for both elements of the contact; E₁, E₂ – the corresponding modulus of elasticity; r₁, r₂ –the corresponding radii of the two spheres.

In the case of a contact between a sphere with the radius r₁ = r and a plane r₂ =$\infty$, as shown in Figure 7b., and the material in both contact elements is the same (σ₁ = σ₂ = σ, E₁ = E₂ = E) then Equation (16) becomes:

$$a=\sqrt[3]{1,5{\left(1-\sigma \right)}^2\frac{F\cdot r}{E}}$$

(17)

or, for r, σ and E assumed constant it can be expressed as:

$$a\approx \sqrt[3]{F}$$

(18)

The dependence between the contact resistance and contact force results:

$$R=\frac{\rho }{2a}\approx \sqrt[3]{F}$$

(19)

Consequently, in the case of plastic deformation, the radius of the contact circle is:

$$a=\sqrt{\frac{F}{\pi H}}$$

(20)

and the dependence is:

$$R=\frac{\rho }{2a}=\frac{\rho }{2}\sqrt{\frac{\pi H}{F}}=\frac{1}{\sqrt{F}}$$

(21)

Technical contacts usually employ plastic deformation. Equation (21) does not show a dependence of resistance on the curvature of the contact, but highlights the fact that the contact resistance is inversely proportional to the contact force at a fractional power. Experience shows that there is a dependence of the contact resistance on the radius of the contact element on silver coated copper contacts.

2.3. The Thermal Regime of Electrical Contacts

From the point of view of the thermal regime of electrical contacts, electric current passing through a contact causes the development of heat (Joule effect) because of contact resistance. In permanent thermal regime, the overtemperature due to contact resistance must be maintained within relatively low limits of 2–10 degrees, whereas in short-term thermal regime such as short-circuit, the overtemperature must not lead to melting of the contact elements. Studying of thermal behavior is carried out for closed contacts, as if the contact elements are separated, additional physical processes such as arcing and material migration occur.

Figure 8 shows a contact point in whose touching area a number of m flow lines (tubes) of current i converge, such that the total current in the touch area is I = m ^. i. Figure 9 schematically details such a tube of current, which in the contact point covers the area ΔA₀.

Calculations are made with the following assumptions:

• − The thermal flux developed in a certain current tube is transmitted outwards only through that current tube. There is no heat transfer between adjacent points m and n in Figure 8, which are assumed to be at the same temperature;
• − The highest overtemperature is at ΔA₀ surface, which also defines an isothermal surface, as a result of which the greatest striction is located in the contact zone;
• − The contact elements are made of the same homogeneous and isotropic material.

By the passing of current i through a current tube, equipotential and isothermal surfaces are defined as shown in Figure 9, considering the surface ΔA₂ with the reference parameters in terms of electric potential V = 0, overtemperature θ = Θ, and absolute temperature T = T_Θ. A surface ΔA_n located far, theoretically at infinity is characterized by V = U/2, where U is the voltage on the contact, θ = 0 and T = T₀.

One can notice that the equipotential surfaces are also isothermal, due to the fact that equal thermal fluxes go through equal thermal resistances. Thus, the equations of electrical resistance R and thermal resistances R_t, for a segment of the current tube between the base surface ΔA₀ and the equipotential surface ΔA₁ situated at a distance d n from ΔA₀ are:

$$\mathrm{d}R=\frac{\mathrm{d}n}{\sigma \cdot \Delta {\overline{A}}_{01}}$$

(22)

$$\mathrm{d}{R}_t=\frac{\mathrm{d}n}{\lambda \cdot \Delta {\overline{A}}_{01}}$$

(23)

where $\Delta {\overline{A}}_{01}$ represents an average value of the area over the distance d n; σ is the electrical conductivity; λ is the thermal conductivity.

For the elementary temperature variation, the equation of Fourier is valid written in the form:

$$\mathrm{d}T=-P_t\cdot \mathrm{d}{R}_t$$

(24)

where P_t is the thermal power and d R_t the thermal resistance traveled by the thermal power. To integrate the differential Equation (24) based on Equations (22) and (23) one can write:

$$\mathrm{d}{R}_t=\mathrm{d}R\frac{\sigma }{\lambda }=\frac{\mathrm{d}R}{\lambda \rho }$$

(25)

which can also be expressed as:

$$\mathrm{d}R=\frac{\mathrm{d}V}{i}$$

(26)

Under these conditions, Equation (25) becomes:

$$-\lambda \rho \mathrm{d}T=V\mathrm{d}V$$

(27)

Equation (27) allows to express the dependence between the absolute temperature of the contact and the electrical potential, of the form:

$$-{\displaystyle {\int }_{T_{\Theta }}^T\lambda \rho \mathrm{d}T}={\displaystyle {\int }_0^{V}V\mathrm{d}V}$$

(28)

To solve Equation (28), λ and ρ are considered as functions of temperature. Two cases can exist, namely:

(a) The overtemperature is moderate, i.e., Θ = 1... 5 degrees. In this case, the product $\lambda \rho$may take an average value $\overline{\lambda \rho }$. Thus, by integrating Equation (28) we obtain:

$$T_{\Theta }-T={\theta }_{\Theta }-\theta =\frac{V^2}{2\overline{\lambda \rho }}$$

(29)

and if the integration is done for a remote area, the contact overtemperature is obtained:

$$T_{\Theta }-T=\Theta =\frac{U^2}{8\overline{\lambda \rho }}$$

(30)

(b) The overtemperature is excessive, meaning that the material tends to melt. In this case, the Wiedemann-Franz-Lorenz law is used:

$$\lambda \rho =LT$$

(31)

By integrating within limits T_Θ, T₀ and, respectively 0, U/2 the Equation (28) it results that:

$$L\left(T_{\Theta }^2-T_0^2\right)=\frac{U^2}{4}$$

(32)

The temperature field of the current paths in which there are contacts is obtained taking into account the permanent end effect, and has the expression:

$$\theta =\left({\theta }^{*}-{\theta }_{\max }\right){\mathrm{e}}^{-\beta x}+{\theta }_{\max }$$

(33)

For the determination of θ^* it is observed that the contact resistance is actually a narrow-place (constriction) resistance and that half of the electrical power developed by the contact resistance is injected into the front surface of the contact, i.e.,:

$$\frac{R{i}^2}{2}=\frac{P}{2}=-\lambda A{\left(\frac{\partial \theta }{\partial x}\right)}_{x=0}$$

(34)

From Equation (33) the temperature gradient at the contact site can be calculated:

$${\left(\frac{\partial \theta }{\partial x}\right)}_{x=0}=-\beta \left({\theta }^{*}-{\theta }_{\max }\right)$$

(35)

so that from (35) the overtemperature of the frontal surface is obtained:

$${\theta }^{*}=\frac{P}{2\lambda A\beta }+{\theta }_{\max }$$

(36)

Note that the equation above has two components, namely θ_c representing the contribution of contact:

$$\theta c=\frac{P}{2\lambda A\beta }=\frac{Ui}{2\lambda A\beta }$$

(37)

and θ_max which is the contribution of the conductor traveled by the current in permanent mode. This means that the contact overtemperature is given by (37) as a function of the power developed by the contact resistance, which implies knowing the voltage U of the contact.

The resistance of the unheated contact (an electric current of zero intensity flows though the contact) differs from the resistance of the heated contact (an electric current of significant intensity flows though it). The dependence between the electrical resistance R_Θ of the heated contact and the electrical resistance R₀ of the unheated contact can be determined using the theory of electric potential. With reference to this, two strictions are considered, namely: striction $S(\lambda ,\rho )$, where $\lambda$ and$\rho$depend on temperature and, striction$S_0({\lambda }_0,{\rho }_0)$, where ${\lambda }_0$ and${\rho }_0$do not depend on temperature. For the two strictions the elementary variation of potential can be expressed as:

$$\mathrm{d}V=\rho \frac{\mathrm{d}n}{\Delta A_n}\mathrm{i};\mathrm{d}{V}_0={\rho }_0\frac{\mathrm{d}n}{\Delta A_n}i$$

(38)

Since the dependence of resistivity on temperature is of the form $\rho ={\rho }_0(1+{\alpha }_R\theta )$, from Equation (38) it can be written that $\mathrm{d}V=\mathrm{d}{V}_0(1+{\alpha }_R\theta )$, and thus:

$$\mathrm{d}{V}_0=\mathrm{d}V\left(1-{\alpha }_R\theta \right)$$

(39)

In this case, Equation (28) can be also written as:

$$-{\displaystyle {\int }_{\Theta }^{\theta }\overline{\lambda \rho }\mathrm{d}\theta }={\displaystyle {\int }_0^{V}V\mathrm{d}V}$$

which becomes, after integration:

$$\theta =\Theta -\frac{1}{2}\cdot \frac{V^2}{\lambda \rho }$$

(40)

where, according to (30),$\Theta =\frac{U^2}{8\overline{\lambda \rho }}$.

Equation (39) with the expression of θ introduced from Equation (40) is integrated within limits Θ, 0 for the heating and 0, U/2 for the potential:

$${\displaystyle {\int }_0^{\frac{U_0}{2}}\mathrm{d}{V}_0}={\displaystyle {\int }_0^{\frac{U}{2}}\mathrm{d}\mathrm{V}}-{\displaystyle {\int }_0^{\frac{U}{2}}\left({\alpha }_R\Theta -\frac{{\alpha }_R{V}^2}{2\overline{\lambda \rho }}\right)\mathrm{d}\mathrm{V}}$$

or:

$$U_0=U-{\alpha }_R\Theta U+\frac{{\alpha }_R{U}^3}{24\overline{\lambda \rho }}=U\left(1-\frac{2}{3}{\alpha }_R\Theta \right)$$

(41)

For ${\alpha }_R\Theta \approx 1$ this can be written as:

$$\frac{U}{U_0}=\frac{1}{1-\frac{2}{3}{\alpha }_R\Theta }\approx 1+\frac{2}{3}{\alpha }_R\Theta$$

(42)

Since the current tube of both strictions is traveled by the same current intensity, it can be written:

$$\frac{U}{U_0}=\frac{R_{\Theta }}{R_0}=1+\frac{2}{3}{\alpha }_R\Theta$$

(43)

and finally:

$$R_{\Theta }=R_0\left(1+\frac{2}{3}{\alpha }_R\Theta \right)$$

(44)

It can therefore be concluded that the resistance of the heated contact can be calculated, starting from the resistance of the unheated contact, with a dependence relationship similar to resistors, in which the temperature coefficient of the resistance is 2/3 of the temperature coefficient of the resistivity.

The Thomson effect refers to the case when an electric current of intensity i travels through a conductor along which there is a temperature gradient, and additional power develops inside the conductor:

$$P_T=\pm k_T\cdot \Delta T\cdot i$$

(45)

where k_T is the Thomson coefficient and ΔT is the temperature difference on the conductor.

The sign ± depends on the direction of current relative to the direction of the temperature gradient. The + sign is considered if the current travels through the conductor in the same direction as the gradient. Current travels through one contact element once in the same direction as the gradient and in the other element in the opposite direction, as shown in Figure 10. As a result, one contact element heats up additionally, and the other cools down. These phenomena manifest themselves in the touch zone, where the temperature is the highest. Accordingly, the anode will be warmer than the cathode.

3. Modeling and Simulation of the Electrical Contacts Using Numerical Methods

We have shown in Section 2 that between the absolute temperature of the contacts and their potential there is a dependence expressed by Equation (28). This equation has a theoretical character, and its integration through analytical methods is difficult. For this reason, the dependence between the potential and temperature of an electrical contact will be analyzed by numerical methods with the aid of the COMSOL Multiphysics software.

3.1. Geometry of the Model

The model of the electrical contact, consisting in an assembly of two contact elements as shown in Figure 11, was built using SolidWorks software. The two elements are identical and the material from which they are built is copper. The geometry and dimensions of one contact element are shown in Figure 12. The advantage of using SolidWorks for the build of the model is given by the dynamic live link feature it has with COMSOL, meaning that any change to the geometry that is done in SolidWorks, is automatically imported into the COMSOL Multiphysics simulation model.

3.2. The Step-by-Step Simulation Setup for Comsol

The simulation was performed for a 3D model in two distinct steps within stationary studies of Structural Mechanics>Solid Mechanics (solid) and Heat Transfer>Electromagnetic Heating>Joule Heating type.

The two contact elements presented above have a cylindrical body with the touching ends shaped similar to a hook. Here the coupling of apparent thermal and electrical resistances by means of mechanical contact pressure takes place. It was considered that the initial temperature of the contact is equal to the ambient temperature, and the difference of electric potential between the contact elements produces heating by Joule effect. Figure 13 shows the geometry of the imported model into COMSOL. In the Definitions menu the Contact Pair type was chosen for the contact elements and their Source Boundaries and Destination Boundaries were defined as shown in Figure 14. Next, for the entire model the material (copper) was chosen with its predefined characteristics from the library as can be seen in Figure 15. The fixed surfaces have been established as shown in Figure 16. These surfaces represent the inner part of the semi-cylindrical areas of the contacts. Since the imported geometric pattern is a cross section through the two electrical contacts, we used the Boundary>Symmetry option from the Physics menu, thus the plane of symmetry of the model has been defined (Figure 17). This approach simplifies computations for half of the model.

The contact method was chosen of Augmented Lagrangian type with the initial value of the contact pressure T_n = 1×10⁷ N/m². Regarding the heat transfer it was picked as Convective heat flux type with a heat transfer coefficient h = 2. The difference of electric potential representing the voltage drop on the simulated contact was imposed by assigning zero voltage (Ground) to the base surface of the cylinder to one of the contact elements, and a voltage of 1 mV (Electric Potential) to the other contact element, as shown in Figure 18a,b.

The finite element mesh is of the Free Tetrahedral type, and the model after meshing looks like shown in Figure 19.

As stated at the start of 3.2 subsection, the simulation is carried out in two stages. In the first stage, only solid mechanics phenomena are simulated using the Solid Mechanics (solids) solver as shown in Figure 20, and in the second stage the Joule effect and heat transfer phenomena are simulated using the Electric Currents and Heat Transfer in Solids solvers as shown in Figure 21. The two stages simulation follows the need to first solve the effect of the contact pressure on the characteristics of the electrical contact elements material, and only next the Joule effect of the contacts and its effect are evaluated.

The simulation of the temperature variation when the current flows through the current path of the contact is done considering the initial temperature of the contact elements T = 293.15 K. The base values of the parameters for the contact elements were setup prior to simulation. Thus, in the Solid Mechanics module of Contact Pair 1, the contact pressure was set to T_n = 1×10⁷ N/m². In the Electric Currents >Pair Electrical Contact module, for the same Contact Pair 1 the Microhardness parameter had an assigned value of H_c = 0.1 GPa. In the Heat Transfer in Solids module, an identical value of the microhardness was set for the Thermal Contact pair.

4. Simulation Results and Discussion

4.1. Results of the Mechanical Simulation

By running the simulation, it is observed that there is an increase in the von Mises stress (Figure 22) in the contact where the contact point and the current path are located. A detailed view of the shape of the von Mises stress can be observed in Figure 23. By analyzing this figure, it is visible that the increase of von Mises stress concentrates to the area of physical contact points between the contact elements.

Figure 24 shows the contact force resulting from the simulation and Figure 25 and Figure 26 shows the detailed variation of this force.

4.2. Results of the Electro-Thermal Simulation

The variation in the electrical potential on the current path of the contact is shown in Figure 27. As expected, it varies between 0 V, which corresponds to the right-side element of the contact which is grounded, and 3 mV on the left side-element. The difference of electric potential is actually the voltage drop on the elements of the electrical contact.

The variation of the electric field amplitude along the current path is shown in Figure 28. The Multislice mode combined with the Streamline Multislice mode was used for this representation. Figure 29 shows a detail of the variation of the electric field amplitude in the touching area of the contact elements.

The variation of temperature on the current path following the simulation is visible in Figure 30, where one can observe that the highest temperature is in the touching area of the electrical contact.

Figure 31 shows the simulation results related to the temperature gradient variation along a cut plane passing through the axis of the current path and parallel to plane XY. It can be seen both from the representation of the magnitude of this quantity <K/m> and from the streamline representation, that the heat transfer from the element of contact with zero electrical potential is greater than that corresponding to the element of contact with an electric potential of 3 mV.

Figure 32 shows the spatial variation of the current path temperature gradient.

Figure 33 shows the variation in the current path temperature by isothermal contours. Isothermal surfaces were drawn for 15 temperature levels. Figure 34 shows the temperature variation of the two contact elements and the variation of the current density on a surface that passes through the axis of symmetry of the current path and is parallel to plane XY, using streamline representation.

A cut plane was built as shown in Figure 35.a., parallel to plane XZ and passing through the origin of the cartesian coordinate system (Y=0). Figures 35.b. and c. show the details regarding the positioning of this cut plane.

In Figure 36 the variation of the current density on the cut plane defined above is presented. In Figure 37 the variation of the current density is shown using a 3-dimensional representation, and the Max/Min-Point diagram showing the density of electric energy. One can notice that the highest density of electric energy is at the touching area of the contact elements.

A dataset element of the type Edge 3D was created (see in Figure 38) along an edge line of one contact element, in order to plot the temperature variation along this edge as presented in Figure 39.

In a similar way, two Edge 3D type datasets are created (Figure 40a,b) for the contours of both electrical contact elements, to be used for the plotting of the variation graph of the electrical potential, as shown in Figure 42.

In order to highlight the variation of the maximum temperature of the current path and its dependence on the value of the contact voltage, the simulations were run for various values of the potential difference between the two elements of the electrical contact. This data along this and the values of the maximum temperature, and the values of the temperature increase as compared to the ambient temperature are all summarized in Table. 1.

This variation of the temperature as a function of the contact voltage values is plotted in Figure 42. As one can observe, the red curve representing this variation has a parabolic shape. A trendline curve of polynomial type 2 represented in black color was also traced with its equation displayed. This trendline overlaps the curve of the variation of the temperature increase.

5. Conclusions

The present paper studies and highlights—using numerical methods of modeling and simulation—the dependence between the contact voltage (the voltage drop on the contact) and the variation of its temperature, but also of the current path within the contact. The study was carried using COMSOL Multiphysics software, specifically three of its simulation modules: solid state mechanics module, electric current module and heat transfer module. The geometric model was build using SolidWorks and imported into COMSOL using the real-time link, which allows for easy modification of the model geometry for other modeling and simulation approaches. The results obtained after the simulations were computed, highlighted the following:

• − The von Misses stress has higher values in the current path area where the flat part intersects with the cylindrical shape for each contact element, and also in the area of the physical touching of the two contact elements;
• − The contact forces have maximum value at the contact points;
• − The electric field amplitude along the current path and the contact elements is maximum at the touching point, where the transverse surface has a minimum value;
• − The temperature gradient along the contact proves that the highest temperature is in the touching area of the electrical contact. The gradient also shows it is higher in the direction from the zero potential contact elements towards the contact point. This supports the theoretical aspects presented in the Introduction section and shown in Figure 10.
• − The variation of temperature in the current path and the contact elements was graphically represented by isothermal surfaces for 15 temperature levels;
• − The dependence between the temperature variation of the two contact elements and the variation of the current density was shown using color legend, on a plane that passes through the axis of symmetry of the current path and is parallel to plane XY, and also using streamline representation;
• − In order to have a more conclusive picture of the current density, we presented its variation on the surface of a horizontal longitudinal cut plane;
• − By drawing the temperature variation diagram along the edges of the contacts, it was observed that its maximum value is found at the touching points of the contact elements takes place;
• − The maximum contact temperature values were computed by simulation for contact voltages between 0 and 5 mV. Thus, the curve of variation of the increase in contact temperature relative to the contact temperature was plotted, finding this variation as parabolic

The study of the obtained values during simulations, indicates an undesirable increase in contact temperature for values higher than 4 mV of voltage drop, with temperatures that may affect the integrity of insulating materials in the electrical contact. This finding is important as damage to the insulation materials may cause the phenomenon of overturning or short circuit in the electrical contact devices, thus leading to fire and damage, aspects that can constitute a future direction of research.