Parallel calculation method of multi-field coupling in three-dimensional electromagnetic relay based on transmission line method

Abstract

The present invention discloses a parallel calculation method of multi-field coupling of a three-dimensional electromagnetic relay based on the transmission line method, and the method is as follows: step one, establishing a three-dimensional electromagnetic-thermal-force multi-field coupling digital prototype model of an aerospace electromagnetic relay; step two, solving the electromagnetic-thermal-force multi-field coupling digital prototype model of the three-dimensional finite element model of the aerospace electromagnetic relay. The present invention mainly solves the dynamic characteristics of the multi-physical field coupling of a three-dimensional aerospace electromagnetic relay in parallel, adopts the edge finite element method to establish the three-dimensional electromagnetic-thermal-force multi-physical field coupling model of the aerospace electromagnetic relay, and introduces TLM into the finite element model solution of the three-dimensional multi-physical field coupling, so as to achieve the purpose of accelerating the solution of dynamic characteristics. The present invention has important theoretical significance for improving the environmental resistance reliability design of aerospace electromagnetic relays, and also has reference significance for shortening the structural design and parameter optimization cycle of aerospace electromagnetic relays.

Description

Technical Field

[0001] The present invention belongs to the field of numerical calculation, and relates to a three-dimensional multi-physical field coupling numerical calculation method suitable for the electromagnetic system of an electrical product, and specifically to a finite element parallel calculation method of three-dimensional multi-physical field coupling based on the Transmission Line Method (TLM). Background Art

[0002] Aerospace electromagnetic relays have the advantages of high conversion depth, excellent physical isolation, and multi-channel synchronous switching. They are widely used in aerospace engineering and are key electronic components that complete functions such as signal transmission, execution control, system power distribution and protection. Their technical level represents the highest level of electromagnetic relays. Compared with other electronic products, aerospace electromagnetic relays are complex electromechanical systems composed of electromagnetic mechanisms and contact units. They often contain magnetic, elastic, and insulating materials, and the working process involves the coupling of multiple physical fields such as electricity, magnetism, heat, force, and fluid. There are many parts and complex processes. Aerospace electromagnetic relays have the worst quality and reliability among all electronic components, especially in terms of environmental resistance and reliability. The differences between different products are obvious, which seriously restricts the development of aerospace engineering.

[0003] Compared with traditional experimental measurements, although the scheme combining digital prototype models with experimental measurements can flexibly and effectively measure and evaluate dynamic characteristics, the existing commercial finite element software has problems such as long time consumption and poor convergence when calculating the dynamic characteristics of three-dimensional digital prototype models, which greatly affects the design and production cycle of aerospace electromagnetic relays. The product design of aerospace electromagnetic relays is divided into three stages: overall design, parameter design, and tolerance design, referred to as three-stage design. Only through three-stage design can the designed products be guaranteed to have the characteristics of long life and high reliability from the source of design. The three-stage design requires continuous adjustment of model parameters. After each adjustment, the digital prototype model needs to be re-solved to optimize the performance parameters. The repeated simulation calculations in the optimization design process take a lot of time. Therefore, the calculation accuracy and solution efficiency of the digital prototype model are crucial for the three-stage design of aerospace electromagnetic relays.

[0004] As the core technology of digital prototype model, finite element method (FEM) is widely used in simulation analysis in the field of electrical appliances. The calculation accuracy and speed of finite element method largely determine the efficiency of electrical appliance optimization design based on digital prototype technology. Compared with two-dimensional problems, the solution method of three-dimensional electromagnetic field is quite different, and the numerical algorithm of three-dimensional multi-physics field coupling model still has the problem of low calculation efficiency. Therefore, it is of great significance to extend the acceleration algorithm to the solution of three-dimensional finite element model. The action process of aerospace electromagnetic relay involves electromagnetic-thermal-force multi-physics field coupling calculation, and the influence of different temperature environments on dynamic characteristics is also very significant. The solution of a single physical field can no longer meet the actual production needs of aerospace electromagnetic relays. However, the digital prototype model of multi-physics field coupling takes more time to solve than the single field model, and its solution efficiency directly affects the design and production cycle of aerospace electromagnetic relays. Therefore, studying the parallel calculation of finite element with multi-physics field coupling has important practical and engineering significance. Summary of the invention

[0005] The purpose of the present invention is to provide a three-dimensional electromagnetic relay multi-field coupling parallel calculation method based on the transmission line method. The method aims to realize the accelerated solution of the digital prototype model of the three-dimensional multi-field coupling of the electrical appliance, and is intended to fill the gap of multidisciplinary cross-cutting and parallel numerical calculation of the electromagnetic-thermal-force multi-physics field coupling modeling of the three-dimensional digital prototype model of the electrical appliance. The present invention mainly solves the dynamic characteristics of the multi-physics field coupling of the three-dimensional aerospace electromagnetic relay in parallel, adopts the edge finite element method to establish the three-dimensional electromagnetic-thermal-force multi-physics field coupling model of the aerospace electromagnetic relay, and introduces TLM into the finite element model solution of the three-dimensional multi-physics field coupling, so as to achieve the purpose of accelerating the solution of the dynamic characteristics. The present invention can be directly applied to the simulation calculation of the three-dimensional digital prototype model of the aerospace electromagnetic relay, which has important theoretical significance for improving the environmental reliability design of the aerospace electromagnetic relay, and also has reference significance for shortening the structural design and parameter optimization cycle of the aerospace electromagnetic relay.

[0006] The objective of the present invention is achieved through the following technical solutions:

[0007] A three-dimensional electromagnetic relay multi-field coupling parallel calculation method based on a transmission line method comprises the following steps:

[0008] Step 1: Establish a three-dimensional electromagnetic-thermal-mechanical multi-field coupling digital prototype model of aerospace electromagnetic relay:

[0009] Step 1: Establish a three-dimensional finite element model of the electromagnetic system of the aerospace electromagnetic relay according to its geometric dimensions, discretize the solution domain into a finite number of tetrahedral units connected in size and shape, and control the number of grids according to the need for solution accuracy to complete the meshing of the finite element model;

[0010] Step 1 and 2: Export and analyze the mesh information of the three-dimensional finite element model to obtain the number of tetrahedral elements, the number of nodes and the area to which the elements belong, as well as the node coordinates and node number mesh information of the tetrahedral elements, and calculate the edge lengths of all tetrahedral elements;

[0011] Step 13: When the coil of the aerospace electromagnetic relay is powered on, a magnetic field will be generated. The transient magnetic field is described by the vector magnetic potential A, and its control equation is:

[0012]

[0013] Where A is the vector magnetic potential; ν is the magnetic resistivity; σ is the electrical conductivity; J e is the coil current density; t is the time;

[0014] Step 14: Use edge finite element method to solve the three-dimensional finite element model:

[0015] Step 141: Split the entire three-dimensional solution region Ω into many smaller tetrahedral units and the edges between tetrahedrons. e The vector magnetic potential vector A of any point in e (x,y,z), approximated by a linear combination of known interpolation functions:

[0016]

[0017] Where i is the number of edges of the tetrahedral unit; W i (x, y, z) is the interpolation function based on the edge vector; is the line integral of the vector magnetic potential field A(x, y, z) along the direction of the i-th edge; the interpolation function W i It is determined only by the four node position coordinates of the tetrahedral element, W i The projection on edge i is 1, and the projection on all other edges is 0;

[0018] Step 142: The edge vector shape function W corresponding to the kth (k = 1, 2, ..., 6)th edge of the tetrahedral element k And its curl is expressed as:

[0019]

[0020] In the formula, l ij represents the edge length between node i and node j, and its positive and negative is determined by the direction vector; N1, N2, N3 and N4 represent the volume coordinates of the tetrahedral unit; ▽N1, ▽N2, ▽N3 and ▽N4 represent the gradient of the volume coordinates of the tetrahedral unit;

[0021] Step 143. The gradient of the volume coordinates of the tetrahedral unit is expressed as:

[0022]

[0023] Where V e is the volume of the tetrahedral unit; is the edge vector of the six edges of the tetrahedral unit;

[0024] Step 144: For a single tetrahedral unit Ω e , where the four nodes are numbered n = 1, 2, 3, 4, the volume V of the tetrahedral unit can be calculated as follows:

[0025]

[0026] In the formula, x n ,y n ,z n They are the x, y, and z coordinates corresponding to the four nodes in the tetrahedron element.

[0027] Step 1, 4, and 5: Use the Galerkin method to compare the residuals with the six weighting functions W i Integrate the product and force the integral to be zero. After considering the Coulomb specification, the unit discrete equation in the form of a 10×10 matrix is ​​obtained:

[0028]

[0029] in

[0030]

[0031]

[0032] Where N i is the volume coordinate of the tetrahedral unit; J e is the coil current density, and the coil current can be obtained through field-circuit coupling; A j is the line integral of the vector magnetic potential field A along the direction of the jth edge;

[0033] Step 146: Establish a three-dimensional transient magnetic field transmission line model to solve the magnetic field equation in step 145. The specific steps are as follows:

[0034] The transmission line method decouples the three-dimensional nonlinear system into a linear network and a nonlinear network, and obtains the equivalent circuit of the transmission line model. The equivalent circuit contains the circuit models of the incident stage and the reflection stage. The two stages are performed alternately in the linear network and the nonlinear network until convergence, thereby completing the solution of the three-dimensional nonlinear system, where:

[0035] (1) The submatrix [K ij ] 6×6 and [D ij ] 6×6 The nonlinear resistance and linear capacitance of the equivalent circuit model are expressed as:

[0036]

[0037] (2) The nonlinear resistor is replaced by the transmission line and the reflection box in the reflection stage. The characteristic conductance of the transmission line is expressed by the following formula:

[0038]

[0039] (3) The linear capacitor is connected in parallel with the nonlinear resistor. The time derivative of the electromagnetic field is related to the capacitor. The capacitor is replaced by the resistor and the current source at the previous moment through the backward Euler method. The calculation formulas for admittance and current source are:

[0040] Y Cij =C ij / Δt,Y Ci =C i / Δt,i,j∈[1,6]

[0041]

[0042]

[0043] (4) In the incident phase, the node voltage A i and branch voltage (A i -A j ) enters the reflection box. In each reflection box, the incident voltage of one port will be scattered to all ports according to the coupling relationship of the nonlinear resistor. After the reflection stage, the incident voltage leaves the reflection box and becomes the reflected voltage. and reflected voltage The nonlinear relationship between is expressed as:

[0044]

[0045] The NR iteration is used to solve this 6×6 system of equations and the Jacobian matrix is ​​calculated according to the chain rule:

[0046]

[0047] After solving the 6×6 equations by NR iteration, the unit reluctance ν is obtained e The value of

[0048] (5) In the incident phase, the reflected voltage containing the unit nonlinear resistance information is returned to the global circuit network, the transmission line and the reflected voltage are replaced by the Norton circuit, and then the node voltage of the global circuit network is updated through the linear network, and the incident voltage of the next iteration is calculated based on the node voltage. In addition, the circuit network in the incident phase is a large-scale linear network assembled by all unit admittances and Coulomb norms in the grid, and remains unchanged during the iteration process of each time step;

[0049] Through steps one, four and six, the vector magnetic potential distribution of the magnetic field can be obtained in parallel;

[0050] Step 147. According to the edge element formula, for each tetrahedral unit, the magnetic induction intensity B is represented by the vector magnetic potential A, and the three-dimensional nonlinear transient magnetic field equation is expressed as:

[0051]

[0052] Step 15: Use the field-circuit indirect coupling scheme to achieve the coupling between the relay electromagnetic field and the circuit:

[0053] Step 151. The equation of the aerospace electromagnetic relay coil is expressed as:

[0054]

[0055] In the formula, u is the coil voltage, i is the coil current, R is the coil resistance, and ψ is the magnetic flux of the excitation coil;

[0056] Step 152: Rewrite the formula in step 151 into the iterative form of a nonlinear equation:

[0057] f (t) =iR+[ψ(t)-ψ(t-1)] / Δt-u=0

[0058] Where ψ(t) is the magnetic flux at the current time step, ψ(t-1) is the magnetic flux at the previous time step, and Δt is the discrete time step;

[0059] Step 153: After obtaining the vector magnetic potential distribution of each edge of the coil unit according to the three-dimensional FEM model, the magnetic flux at the current time step is calculated by the following formula:

[0060]

[0061] Among them, N coil is the number of coil turns, S coil is the cross-sectional area of ​​the coil, n is the unit vector along the direction of the coil current, V coil is the volume of the coil area unit, N is the number of edges;

[0062] Step 154: Use Newton's iteration method to construct the formula in step 152 into the following iterative format to calculate the current value i at the (t+1) moment: (t+1) , the iteration format is as follows:

[0063]

[0064] Get the current value i at the next moment (t+1) By inputting the transient magnetic field for iterative calculation, the vector magnetic potential value at the next moment is obtained, and the field-circuit coupling between the electromagnetic field of the aerospace electromagnetic relay and the control circuit can be realized;

[0065] Step 16: Use the three-dimensional edge element and the unit load method to derive the parallel calculation method of electromagnetic force and electromagnetic torque:

[0066] Step 161. According to the unit load method theory, the electromagnetic force is expressed as:

[0067]

[0068] Where q is the displacement direction of the armature, F q W is the electromagnetic force on the armature in the q direction. mag is the magnetic field energy;

[0069] Step 162: The electromagnetic force at node k in the x direction in tetrahedral element e is described as:

[0070]

[0071] Step 163: Combined with the three-dimensional nonlinear transient magnetic field equation, the unit node force F ek Simplified to:

[0072]

[0073] In the formula, It can be obtained through the BH curve of ferromagnetic materials. It can be obtained by the following formula:

[0074]

[0075] Step 164. Volume V of tetrahedral unit e x k The partial derivative of is expressed as:

[0076]

[0077] Step 165: After calculating the electromagnetic force in the x direction, calculate the electromagnetic forces in the y and z directions in the same way; the electromagnetic force and electromagnetic torque are calculated in parallel;

[0078] Step 17: Solve the electromagnetic field-mechanical motion coupling problem:

[0079] Step 171. After obtaining the electromagnetic force of each node, calculate the electromagnetic torque of each node according to the following formula:

[0080]

[0081] Where M e is the electromagnetic torque of the e-th unit, is the magnitude of the electromagnetic force at the node, is the vertical distance from the action point to the rotation axis, N n is the number of nodes in each unit;

[0082] Step 172: Accumulate the electromagnetic torque of all contour nodes to obtain the total electromagnetic torque:

[0083]

[0084] Where M mag is the electromagnetic torque acting on the armature, N e is the number of tetrahedral units;

[0085] Step 173. After calculating the electromagnetic torque, the dynamic equation of the armature of the aerospace electromagnetic relay in the rotational motion is described as:

[0086]

[0087] Where I is the moment of inertia, α is the angular acceleration of the armature, θ is the rotation angle of the armature, and M r is the mechanical resistance torque acting on the armature;

[0088] Step 18: Finite element calculation of thermal conduction of aerospace electromagnetic relay:

[0089] Step 181. The three-dimensional partial differential equation of the transient temperature field is expressed as:

[0090]

[0091] In the formula, ρ is the density of the material, C is the specific heat capacity of the material, T is the temperature of the solid or air, λ is the thermal conductivity of the medium, and Q is the heat generation power per unit volume;

[0092] Step 182: Use the tetrahedron node temperature distribution to characterize the temperature field, and use the FEM discretization process to obtain the following equation:

[0093]

[0094] After the above formula is discretized in time using the backward Euler method, the following algebraic equation is obtained:

[0095]

[0096] In the formula,

[0097]

[0098] Where l,n=1,2,3,4; is the density of the material in the unit area, C e is the specific heat capacity of the material in the unit area, λ e is the thermal conductivity of the unit medium, G e is the internal heat source intensity, p l ,q l 、r l and l is a function related to the node coordinates x, y and z;

[0099] Step 183. The heat sources in aerospace electromagnetic relays mainly include Joule loss and eddy current loss, of which eddy current loss P eddy After the magnetic field calculation is completed, it is obtained through post-processing and expressed as:

[0100]

[0101] Where E is the electric field strength, Ω core It is the ferromagnetic region;

[0102] Step 2: Solve the electromagnetic-thermal-mechanical multi-field coupled digital prototype model of the three-dimensional finite element model of the aerospace electromagnetic relay:

[0103] Step 21: Read the temperature distribution at the last moment and update the material properties of each area;

[0104] Step 22: Read the vector magnetic potential distribution and magnetic reluctance at the current moment, and calculate the Joule loss and eddy current loss;

[0105] Step 2: Calculate the temperature field to update the temperature distribution, and update the conductivity of the coil and the ferromagnetic material;

[0106] Step 24: Calculate the electromagnetic force and electromagnetic torque on the armature according to the unit load method;

[0107] Step 25: Calculate the mechanical motion equation to determine whether the armature rotates. If so, update the mesh information of the three-dimensional electromagnetic-thermal-mechanical multi-field coupled digital prototype model. If the armature remains stationary, there is no need to update the mesh information.

[0108] Step 26: Calculate the coil current value of the next time step based on the field-circuit coupling, and save the vector magnetic potential and magnetic reluctance of the next time step.

[0109] Compared with the prior art, the present invention has the following advantages:

[0110] (1) Two-dimensional models can only be applied to symmetrical structures or specific occasions, but some complex three-dimensional geometric structures cannot be represented in two-dimensional space. Therefore, the present invention uses the edge finite element method to establish a three-dimensional electromagnetic-thermal-mechanical multi-physics field coupling model of aerospace electromagnetic relays. Compared with the single-field digital prototype model, the multi-physics field coupled three-dimensional digital prototype model can comprehensively consider the influence of ambient temperature and internal temperature rise on the dynamic characteristics of aerospace electromagnetic relays, and is closer to the actual working conditions of aerospace electromagnetic relays. At the same time, the method of combining three-dimensional digital modeling technology with experimental measurement provides a new and more economical and efficient solution for testing and evaluating the dynamic characteristics of electromagnetic relays.

[0111] (2) The three-dimensional digital prototype model proposed in this invention realizes the electromagnetic-thermal-mechanical multi-physics field coupling of aerospace electromagnetic relays, provides a detailed implementation plan for the three-dimensional multi-physics field coupling and parallel solution of aerospace electromagnetic relays, and analyzes the influence of temperature rise on the dynamic characteristics of electromagnetic relays. The dynamic characteristics evaluation results of the three-dimensional digital prototype model are very consistent with the experimental measurement results, and the relative error is controlled within 4%. The method combining the three-dimensional digital prototype model and experimental measurement provides a more flexible and efficient solution for measuring and evaluating the dynamic characteristics of aerospace electromagnetic relays.

[0112] (3) When solving a 3D multi-physics field coupled finite element model, the time required to calculate the dynamic characteristics using commercial finite element software is long due to the large number of grids in the 3D model. The present invention extends the TLM parallel algorithm to the solution of 3D multi-physics field coupling problems, decouples the nonlinear hybrid network into a linear network and a nonlinear network, and realizes large-scale parallel computing. In comparison, TLM can improve the calculation speed by up to 2.71 times when solving the dynamic characteristics of three-dimensional multi-physics field coupling, which proves the correctness and effectiveness of the proposed three-dimensional multi-physics field coupling dynamic characteristics acceleration method, which will be very helpful for the rapid measurement and evaluation of the dynamic characteristics of aerospace electrical appliances. TLM is very suitable for solving three-dimensional large-scale nonlinear problems. With the improvement of the parallel computing power of computer hardware, the efficiency of TLM in solving three-dimensional multi-physics field coupling problems can be further improved, and it has good application prospects.

[0113] (4) The three-dimensional digital prototype modeling method of the present invention can be applied to relays and contactors of other different structures, and can be widely used in engineering practices such as optimization design and electrical life test of electrical products, with good versatility and practicality. The proposed parallel computing method can shorten the evaluation time of the dynamic characteristics of electrical appliances, thereby shortening the production and optimization cycle of products and reducing time costs. In addition, the present invention can transform achievements and develop professional three-dimensional finite element simulation software in the field of electrical appliances, breaking the monopoly of foreign countries on finite element simulation software, and further highlighting the beneficial effects and application value of the present invention. BRIEF DESCRIPTION OF THE DRAWINGS

[0114] Figure 1 is a first-order tetrahedral edge element;

[0115] Figure 2 Control circuit for aerospace electromagnetic relays;

[0116] Figure 3 Transmission line equivalent circuit of three-dimensional nonlinear system, (a) equivalent circuit of nonlinear unit equation, (b) transmission line model is replaced by Norton circuit in the incident stage, (c) nonlinear resistor is replaced by coupled nonlinear unit model (green box) and transmission line in the reflection stage;

[0117] Figure 4 A massively parallel implementation of a detailed 3D-TLM-FEM;

[0118] Figure 5 This is the numerical calculation flow chart of electromagnetic-thermal-mechanical multi-field coupling for the three-dimensional finite element model of aerospace electromagnetic relay;

[0119] Figure 6 Pseudocode for numerical calculation of electromagnetic-thermal-mechanical multi-field coupling for three-dimensional finite element model of aerospace electromagnetic relay;

[0120] Figure 7 This is a physical picture of aerospace electromagnetic relay;

[0121] Figure 8 It is the BH curve of the nonlinear material DT4C in the ferromagnetic region;

[0122] Fig. 9 The magnetic flux density distribution diagram of aerospace electromagnetic relay with specified sampling time, (a) t1 = 1.5 ms, (b) t2 = 4.5 ms;

[0123] Fig.10 It is a high-speed video test platform for dynamic characteristics of aerospace electromagnetic relays;

[0124] Fig.11 This is the dynamic characteristic curve of aerospace electromagnetic relay;

[0125] Fig.12 It is a temperature rise test platform for aerospace electromagnetic relays;

[0126] Fig.13 The change of coil current and internal temperature of aerospace electromagnetic relay over time;

[0127] Fig.14 The temperature distribution of the internal components of the aerospace electromagnetic relay after repeated operation for 20 minutes;

[0128] Fig.15 The dynamic characteristics of the aerospace electromagnetic relay after repeated action for 20 minutes;

[0129] Fig.16 Solve information for TLM dynamics using different meshes and thread counts. DETAILED DESCRIPTION

[0130] The technical solution of the present invention is further described below in conjunction with the accompanying drawings, but is not limited thereto. Any modification or equivalent replacement of the technical solution of the present invention without departing from the spirit and scope of the technical solution of the present invention should be included in the protection scope of the present invention.

[0131] The present invention provides a three-dimensional electromagnetic relay multi-field coupling parallel calculation method based on the transmission line method. First, a three-dimensional finite element model of the electromagnetic system is established according to the geometric structure and parameters of the aerospace electromagnetic relay; then, the differential equations of the electromagnetic field, circuit, electromagnetic torque, mechanical motion and temperature field are established according to the physical characteristics of the aerospace electromagnetic relay, and discretized according to the Galerkin finite element method; finally, appropriate intermediate variables are selected to realize the indirect coupling of the electromagnetic field, circuit, electromagnetic force, mechanical motion and temperature field equations, and complete the establishment of a three-dimensional electromagnetic-thermal-force multi-physical field coupling digital prototype model. TLM decouples the nonlinear transient magnetic field into a linear network and a nonlinear network. The transmission line admittance matrix of the linear network remains unchanged during the iteration process, and the nonlinear network can be equivalent to multiple independent Norton circuits for parallel calculation, thereby realizing the accelerated solution of the dynamic characteristics of the aerospace electromagnetic relay. The flowchart and pseudo code of the parallel numerical calculation are respectively as follows. Figure 5 and Figure 6 As shown, the specific execution steps are as follows:

[0132] 1. Establish a three-dimensional electromagnetic-thermal-mechanical multi-field coupling digital prototype model of aerospace electromagnetic relay. The specific modeling process is as follows:

[0133] Step 1: Using the commercial finite element software COMSOL The software's three-dimensional mapping function establishes a three-dimensional finite element model based on the geometric dimensions of the electromagnetic system of the aerospace electromagnetic relay, discretizes the solution domain into a finite number of connected tetrahedral units of finite size and shape, and controls the number of grids according to the need for solution accuracy to complete the meshing of the finite element model. The mesh information of the three-dimensional finite element model is exported and parsed to obtain the mesh information including the number of tetrahedral units, the number of nodes, the area to which the units belong, the node coordinates and node numbers of the tetrahedral units, and calculate the edge lengths of all tetrahedral units.

[0134] Step 2: When the coil of the aerospace electromagnetic relay is powered on, a magnetic field will be generated. The transient magnetic field can be described by the vector magnetic potential A, and its control equation is:

[0135]

[0136] Where A is the vector magnetic potential; ν is the magnetic resistivity; σ is the electrical conductivity; J e is the coil current density; t is the time.

[0137] The traditional finite element method defines the unknown function at the unit node. When the dual curl equation is discretized using the vector magnetic potential method according to the finite element method, since each three-dimensional finite element model has three components x, y, and z at the node, the number of finite element equations finally obtained is three times the number of nodes, which will bring a great burden to the computer memory. The edge finite element method uses the product of the tangential component of the field vector along the edge and the edge length as the unknown function. This ensures that the tangential component of the magnetic field intensity H of adjacent units sharing the same edge is equal, does not constrain the normal component, and avoids the problem of large memory consumption caused by the use of the three-dimensional vector magnetic potential method. Therefore, the present invention adopts the edge method to solve the three-dimensional finite element model. The finite element method splits the entire three-dimensional solution area Ω into many smaller tetrahedral units and edges between tetrahedrons. For example Figure 1 As shown, for the tetrahedral unit Ω e The vector magnetic potential vector A of any point in e (x,y,z), which can be approximated by a linear combination of known interpolation functions

[0138]

[0139] Where i is the number of edges of the tetrahedral unit; W i (x, y, z) is the interpolation function based on the edge vector; is the line integral of the vector magnetic potential field A(x, y, z) along the direction of the i-th edge. Interpolation function W i It is determined only by the four node position coordinates of the tetrahedral element, W i The projection on edge i is 1, and the projection on all other edges is 0.

[0140] The order of the six edges and node directions of the tetrahedral element are as follows Figure 1 As shown, the edge vector shape function W corresponding to the kth (k = 1, 2, ..., 6) edge k and its curl can be expressed as:

[0141]

[0142] In the formula, l ij represents the edge length between node i and node j, and its positive and negative is determined by the direction vector; N1, N2, N3 and N4 represent the volume coordinates of the tetrahedral unit; ▽N1, ▽N2, ▽N3 and ▽N4 represent the gradient of the volume coordinates of the tetrahedral unit.

[0143] The gradient of the volume coordinates of the tetrahedral unit can be expressed as:

[0144]

[0145] Where V e is the volume of the tetrahedral unit; and for Figure 1 The edge vectors of the six edges of the tetrahedral element shown.

[0146] For a single tetrahedral unit Ω e , whose four nodes are numbered n = 1, 2, 3, 4, then the volume V of the tetrahedral unit is e It can be calculated by the following formula:

[0147]

[0148] In the formula, x n ,y n ,z n They are the x, y, and z coordinates corresponding to the four nodes in the tetrahedron element.

[0149] The Galerkin method is a weighted residual method that combines the residual with six weighting functions W i Integrate the product of and force the integral to zero. After considering the Coulomb norm, the unit discrete equation in the form of a 10×10 matrix is ​​obtained:

[0150]

[0151] in

[0152]

[0153]

[0154]

[0155]

[0156]

[0157] Where N i is the volume coordinate of the tetrahedral unit; J e is the coil current density, and the coil current can be obtained through field-circuit coupling; A j is the line integral of the vector magnetic potential field A along the jth edge; variable λ e It has no special physical meaning.

[0158] In addition, the magnetic reluctance ν e It is related to the magnetic induction intensity B, and B is related to A i This means that equation (6) is A i The nonlinear equations are generally assembled and the Newton-Raphson (NR) iteration method is often used to solve the nonlinear equations. According to the edge element formula, for each tetrahedral unit, the magnetic induction intensity B can be expressed by the vector magnetic potential A:

[0159]

[0160] Wherein, the values ​​of i1 and i2 are shown in formula (3).

[0161] However, the above discussion only gives six local A i In the entire solution domain, the degrees of freedom of the six edges of all tetrahedral elements must be constrained to each other and subject to Maxwell's equations. Each edge has only one global direction and value A. k When an edge is shared by multiple tetrahedral elements, all local A i Should be equal to A k . All local W i The signs of the components on the edges are adjusted so that they all point in the global direction, which ensures the continuity of the tangent A. Therefore, the normal B is continuous between elements.

[0162] During the assembly process, a global matrix is ​​formed, and each edge has only one global unknown (the multiple unit unknowns that share the edge correspond to the global unknown). The unit matrix and vectors are directly added to the global matrix or vector at the corresponding position. Therefore, all constraints between the local unit and the global matrix are satisfied, and after adding appropriate boundary conditions, the global matrix can be solved.

[0163] However, for three-dimensional problems, the global sparse matrix may have hundreds of thousands of unknowns (the total number of edges in the solution domain). In addition, when the NR iteration method is used for solving, the global matrix needs to be updated and assembled in each iteration, which will lead to excessive memory usage and inefficient solution efficiency. Therefore, the present invention adopts a transmission line parallel algorithm to solve the nonlinear transient magnetic field to improve the solution efficiency of the three-dimensional finite element model.

[0164] Step 3: When the aerospace electromagnetic relay is working, it needs to be powered by an external power supply to generate a magnetic field. The present invention adopts a field-circuit indirect coupling scheme to achieve the coupling between the relay electromagnetic field and the circuit. Figure 2 As shown, the equation of the coil circuit can be expressed as:

[0165]

[0166] Where u is the coil voltage, i is the coil current, R is the coil resistance, and ψ is the magnetic flux of the excitation coil.

[0167] Formula (13) can be rewritten as the iterative form of nonlinear equation:

[0168]

[0169] In the formula, ψ(t) is the flux linkage of the current time step, ψ(t-1) is the flux linkage of the previous time step, and Δt is the discrete time step. After obtaining the vector magnetic potential distribution of each edge of the coil unit according to the three-dimensional FEM model, the flux linkage of the current time step can be calculated by the following formula:

[0170]

[0171] Among them, N coil is the number of coil turns, S coil is the cross-sectional area of ​​the coil, n is the unit vector along the direction of the coil current, V coil is the volume of the coil area unit, and N is the number of edges.

[0172] The Newton iteration method can be used to construct equation (14) into the iterative format of equation (16) to calculate the current value i at the (t+1) moment: (t+1) , the iteration format is as follows:

[0173]

[0174] Get the current value i at the next moment (t+1) The transient magnetic field can then be substituted into the iterative calculation to obtain the vector magnetic potential value at the next moment, thus realizing the field-circuit coupling between the electromagnetic field of the aerospace electromagnetic relay and the control circuit.

[0175] Step 4: After the coil of the aerospace electromagnetic relay is energized, the electromagnetic system generates an electromagnetic torque. When the electromagnetic torque is greater than the counter torque that maintains the contact open, the contact can be closed. The present invention uses a combination of three-dimensional edge element and unit load method to derive a parallel calculation method for electromagnetic force and electromagnetic torque.

[0176] According to the unit load method theory, the electromagnetic force can be expressed as:

[0177]

[0178] Where q is the displacement direction of the armature, F q W is the electromagnetic force on the armature in the q direction. mag is the magnetic field energy, and its calculation formula can be expressed as:

[0179]

[0180] Where W e is the unit magnetic field energy, N e is the number of tetrahedral elements that include nodal forces.

[0181] The electromagnetic force at node k in the x-direction in the tetrahedral element e can be described as:

[0182]

[0183] Combined with formula (12), the unit node force F ek It can be simplified to:

[0184]

[0185] In the formula, can be obtained from the BH curve of ferromagnetic materials, and It can be obtained by the following formula:

[0186]

[0187] The volume V of the tetrahedral element e x k The partial derivative of can be expressed as:

[0188]

[0189] After calculating the electromagnetic force in the x direction, the electromagnetic force in the y and z directions can be calculated in the same way. Based on the above unit load method theory, after calculating the transient electromagnetic field to obtain the vector magnetic potential A and the unit magnetic induction intensity B, the armature electromagnetic force can be calculated. In addition, since the calculations between units are independent of each other, the electromagnetic forces in different units can be calculated in parallel and accumulated in parallel to obtain the electromagnetic resultant force.

[0190] Step 5: The electromagnetic mechanism of the aerospace electromagnetic relay contains mechanical moving parts, and the electromagnetic force will cause the moving parts to move. Therefore, it is necessary to solve the electromagnetic field-mechanical motion coupling problem.

[0191] After obtaining the electromagnetic force of each node, the electromagnetic torque of each node can be calculated according to the following formula:

[0192]

[0193] Where M e is the electromagnetic torque of the e-th unit, F i e is the magnitude of the electromagnetic force at the node, is the vertical distance from the action point to the rotation axis, N n is the number of nodes for each element.

[0194] By adding up the electromagnetic moments of all contour nodes, the total electromagnetic moment can be obtained:

[0195]

[0196] Where M mag is the electromagnetic torque acting on the armature, N e is the number of tetrahedral units.

[0197] After calculating the electromagnetic torque, the dynamic equation of the armature of the aerospace electromagnetic relay in rotational motion can be described as:

[0198]

[0199] Where I is the moment of inertia, α is the angular acceleration of the armature, θ is the rotation angle of the armature, and M r is the mechanical resistance torque acting on the armature.

[0200] It is worth noting that when the armature rotates, the surrounding air domain will deform and the re-meshing technology needs to be performed. However, after the re-meshing is performed, the mesh information before and after the armature movement will change, and the number of units and nodes will also change to a certain extent. Therefore, after the re-meshing is performed, the physical information on the mesh before and after the movement needs to be transferred by mesh mapping to calculate the time-related partial derivatives.

[0201] Step 6: The coil and core areas of the aerospace electromagnetic relay will generate Joule loss and eddy current loss respectively. These losses will be conducted inside the relay in the form of heat energy, causing the temperature of the internal components of the relay to rise. The influence of temperature should not be ignored in the three-stage design process, especially for aerospace standards, which require electrical products to work reliably in low and high temperature environments. Aerospace electromagnetic relays usually adopt a sealed structure to prevent the external environment from affecting its internal components, and its outer shell has good mechanical strength and corrosion resistance. Because the sealing structure restricts the flow of air, heat is difficult to dissipate, especially when the aerospace electromagnetic relay is running for a long time, the internal temperature rise is very obvious. The main form of heat exchange inside the aerospace electromagnetic relay with a sealed structure is heat conduction, and the influence of heat convection and heat radiation on the dynamic characteristics of the electrical appliance is small and can be ignored. Therefore, the present invention mainly studies the finite element calculation of thermal conduction of aerospace electromagnetic relays.

[0202] Heat energy will be conducted through different media inside the electromagnetic relay, then transferred to the surface of the electromagnetic relay device, and finally dissipated through heat exchange with the surrounding medium. Fourier's law is an important law for analyzing heat conduction phenomena. The three-dimensional partial differential equation of the transient temperature field can be expressed as:

[0203]

[0204] In the formula, ρ is the density of the material, C is the specific heat capacity of the material, T is the temperature of the solid or air, λ is the thermal conductivity of the medium, and Q is the heat generation power per unit volume. Material parameters can be found in the material library of commercial finite element software.

[0205] Unlike the vector magnetic potential of the three-dimensional electromagnetic field, temperature is a scalar. The temperature field can be characterized by the node temperature distribution of the tetrahedron. The following equation can be obtained by using the discretization process of FEM:

[0206]

[0207] In the formula,

[0208]

[0209] Where l,n=1,2,3,4; is the density of the material in the unit area, C e is the specific heat capacity of the material in the unit area, λ e is the thermal conductivity of the unit medium, G e is the internal heat source intensity. l ,q l 、r l and l It is a function related to the node coordinates x, y and z, and its functional relationship is as follows:

[0210]

[0211] After time discretization of equation (27) using the backward Euler method, the following algebraic equation can be obtained:

[0212]

[0213] It is worth noting that the boundary conditions of the temperature field are ubiquitous and cannot be ignored in the three-dimensional temperature field analysis. The present invention adopts the first type of boundary conditions, that is, the temperature value on the boundary of the three-dimensional structure is known.

[0214] The dynamic characteristics of aerospace electromagnetic relays will produce a changing magnetic field, and the conductors such as the armature and yoke inside the relay will produce eddy current losses in the changing magnetic field. Eddy current losses will increase the temperature of the core area in the form of heat energy, thus forming local hot spots. It is worth noting that although Joule heat dominates heat conduction, since aerospace electrical products need to be repeatedly opened and closed hundreds of thousands of times during electrical life tests, the eddy current phenomenon caused by repeated magnetic field changes on the temperature rise of aerospace electromagnetic relays should not be ignored, which may affect the electrical life and reliability of aerospace electromagnetic relay products.

[0215] The heat sources in aerospace electromagnetic relays mainly include Joule loss and eddy current loss, of which eddy current loss P eddy After the magnetic field calculation is completed, it is obtained through post-processing and can be expressed as:

[0216]

[0217] Where E is the electric field strength, Ω core It is a ferromagnetic region.

[0218] Based on the above multi-physics field theories and methods, the electromagnetic-thermal-mechanical multi-physics field coupling of the three-dimensional finite element model of the electromagnetic mechanism of the aerospace electromagnetic relay can be realized by indirect coupling.

[0219] 2. Three-dimensional nonlinear finite element parallel calculation based on transmission line method.

[0220] For nonlinear problems, the unit equations are generally assembled globally and then NR iteration is used to solve the nonlinear equations. However, NR iteration requires repeated updates of the Jacobian matrix during each iteration, and the matrix dimension of three-dimensional problems is huge, which makes its solution inefficient. The transmission line method can decouple the hybrid network into linear and nonlinear networks and solve them separately. The solution of the linear network does not require repeated updates of the Jacobian matrix, and the nonlinear network can achieve large-scale parallel calculations, thereby realizing the parallel calculation of the three-dimensional multi-physics field coupled finite element model of the aerospace electromagnetic relay.

[0221] The TLM iteration consists of an incident phase and a reflection phase. Figure 3 (a) shows the submatrix [K ij ] 6×6 and [D ij ] 6×6 The nonlinear resistance and linear capacitance of the equivalent circuit model. The six red nodes represent the circuit nodes, and the voltage value A on the node i It is equal to the projection value of the vector magnetic potential A in the direction of the corresponding edge i. The blue line and the yellow line represent the ground conductance and mutual conductance at the six circuit nodes respectively. The circuit elements of each branch are shown on the right, and their values ​​are defined as follows:

[0222]

[0223]

[0224]

[0225] In the formula, G i is the conductivity of edge i to ground, G ij is the mutual conductance between edges i and j. e By A i Determine the unknown reluctance, ν e The value of C will be obtained through NR iteration in the reflection stage. ij and C i is the equivalent capacitance, which is related to the time derivative term.

[0226] exist Figure 3 In the reflection stage shown in (c), the nonlinear resistor is replaced by a transmission line and a reflection box. The characteristic conductance of the transmission line can be expressed as follows:

[0227]

[0228]

[0229] In the formula, is the guessed value of the unit magnetic reluctance. The selection of this value is arbitrary, but in order to make TLM converge quickly, the selection of this value should be as close to the true value v as possible. e .

[0230] like Figure 3 As shown in (a), the linear capacitor is connected in parallel with the nonlinear resistor. The time derivative of the electromagnetic field is related to the capacitor, and the capacitor can be replaced by the resistor and the current source at the previous moment through the backward Euler method. The calculation formula for admittance and current source is:

[0231]

[0232]

[0233]

[0234] In the formula, the superscript (t) represents the current time step, and (t-1) represents the previous time step.

[0235] like Figure 3 As shown in (c), the node voltage A i and branch voltage (A i -A j ) enters the reflection box. In each reflection box, the incident voltage at one port is scattered to all ports according to the coupling relationship of the nonlinear resistor. After the reflection stage, the incident voltage leaves the reflection box and becomes the reflected voltage. For example, Entering port 1, and distributed and guided to ports 1-21 according to the relationship between G1 and other conductances, all incident voltages are recombined to form reflected voltages. e The node voltage A i Rather than the node voltage difference (A i -A j ), so only the blue ports participate in the NR iteration. and reflected voltage The nonlinear relationship between can be expressed as:

[0236]

[0237] Among them, Y Gi is a fixed guess value, but G i yes The function is determined by equation (12) and the BH curve. The NR iteration can be used to solve this 6×6 system of equations and calculate the Jacobian matrix according to the chain rule:

[0238]

[0239] After solving the 6×6 equations through NR iteration, the unit reluctance ν can be obtained e The value of ν e G for calculating ports 7-21 ij The reflected voltage of port 7-21 can be obtained by the following formula:

[0240]

[0241] The above NR iteration solution process does not involve the node voltage difference (A i -A j ), the entire iteration process only involves the node voltage A i , so there is no need to perform additional transformation and rewriting on formula (12).

[0242] In such Figure 3 In the incident phase shown in (b), the reflected voltage containing the unit nonlinear resistance information is returned to the global circuit network. The transmission line and the reflected voltage are replaced by the Norton circuit, and then the node voltage of the global circuit network is updated through the linear network, and the incident voltage of the next iteration is calculated based on the node voltage. In addition, the circuit network in the incident phase is a large-scale linear network assembled from all the unit admittances and Coulomb norms in the grid, and it remains unchanged during the iteration process of each time step.

[0243] In summary, in the reflection stage, the incident voltage based on the edge variables will be injected into each tetrahedral unit. Since the NR iterations of each tetrahedral unit are independent of each other and the scale of the NR iteration is very small (6×6), large-scale parallel calculation of the reflected voltage can be achieved. In addition, TLM only needs to calculate the global admittance matrix once at each time step, and the admittance matrix will not be changed during the iteration process. The traditional NR iteration method needs to update the Jacobian matrix during each nonlinear iteration, resulting in large amount of calculation and low computational efficiency. Therefore, TLM parallel computing can have a good acceleration effect in 3D problems with a large number of units. The detailed large-scale parallel implementation plan based on 3D-TLM-FEM is shown below. Figure 4 shown.

[0244] By applying TLM to the transient electromagnetic field solution of the three-dimensional FEM model and adopting a parallel scheme in the calculation process of the electromagnetic force and electromagnetic torque, the multi-physics field coupling parallel calculation of the aerospace electromagnetic relay can be realized.

[0245] Example:

[0246] 1. Example description and parameter setting

[0247] The actual picture and structure of the aerospace electromagnetic relay are as follows Figure 7 As shown. According to the structural diagram of the electromagnetic system of the aerospace electromagnetic relay, its finite element model is established. The three-dimensional finite element model of this type of aerospace electromagnetic relay contains 24,337 nodes, 141,560 units and 167,107 edges at the initial position of the armature. The movement mode of this type of aerospace electromagnetic relay is rotational. During the rotation of the armature, the geometric information of the FEM model will be updated and the re-meshing will be performed. The three-dimensional finite element model at the armature closure position contains 27,843 nodes, 162,550 units and 191,613 edges. In addition, the dynamic characteristics measurement and temperature rise test are carried out at room temperature 298.15K. The rated voltage of this type of aerospace electromagnetic relay is 27V, the number of coil turns is 4950 turns, and the coil resistance is 560Ω.

[0248] The dynamic characteristics of aerospace electromagnetic relays are nonlinear because the material properties of the internal components of the relays change nonlinearly. For example, the nonlinear BH curve of the ferromagnetic region in the electromagnetic mechanism is as follows: Figure 8 shown.

[0249] TLM and COMSOL The calculations are all performed on a 64-bit Windows server with 64GB of memory, including an AMD Ryzen 97950X3D CPU (4.20GHz, 16 cores and 32 threads) and an NVIDIA GeForce GT 1030 (4GB) GPU. In order to give full play to the performance advantages of multi-core CPUs and improve the efficiency of TLM parallel algorithms in solving multi-physics field coupled finite element models, the present invention uses multi-core CPUs for parallel computing, and parallelizes the TLM admittance calculations of tetrahedral units, the incident voltage and reflected voltage calculations, and the post-processing calculations of electromagnetic forces, magnetic flux, etc. Although pre-memory allocation and data transfer will take up a certain amount of time and affect the computing efficiency, the three-dimensional finite element model has a large number of grids and a long time is required for single-step solution. Therefore, multi-core parallel computing is more advantageous than serial computing. At the same time, in the commercial finite element software COMSOL The dynamic characteristics simulation of aerospace electromagnetic relays was also carried out to compare with the numerical method proposed in this invention. In addition, for the numerical calculation and simulation of dynamic characteristics, the maximum angle allowed for the armature to rotate is 7.2°, the total simulation time is set to 5ms, and the time step is set to 0.1ms. For the numerical calculation and simulation of the temperature rise test, the total simulation time is set to 20 minutes, and the time step is set to 1s.

[0250] 2. Analysis of the results of the embodiment

[0251] (1) Measurement and evaluation of dynamic characteristics

[0252] The dynamic characteristics of the electromagnetic mechanism of aerospace electromagnetic relays mainly include the relationship between coil current, coil flux, armature rotation angle, angular velocity and angular acceleration, electromagnetic torque, etc. and time. After the coil is energized, the coil current will gradually increase, and the armature will rotate under the drive of the electromagnetic torque. Under conditions of long-term power-on or repeated power-on, the calculation of the transient electromagnetic field needs to be coupled with temperature changes. The temperature rise of the internal components of the electromagnetic relay will affect the changes in material properties and affect the calculation of dynamic characteristics. Therefore, the dynamic characteristics of aerospace electromagnetic relays involve electromagnetic-thermal-mechanical multi-physical field coupling calculations. The present invention uses TLM to solve the dynamic characteristics of aerospace electromagnetic relays in parallel. Fig. 9 The magnetic induction intensity distribution cloud diagrams at sampling times t1=1.5ms and t2=4.5ms are given.

[0253] from Fig. 9 It can be seen that at time t1 = 1.5ms, the armature made a preliminary movement, the coil current was still in the rising stage, and the maximum value of the magnetic induction intensity was 1.73T. As time went on, the armature rotated and finally reached the pull-in position. At time t2 = 4.5ms, the coil current reached the steady-state value and the magnetic field tended to be stable. The magnetic induction intensity reached saturation and showed a highly nonlinear distribution, with a maximum value of 2.04T.

[0254] In order to further verify the accuracy of the three-dimensional digital prototype model proposed in this invention in calculating the dynamic characteristics of aerospace electromagnetic relays, COMSOL The dynamic characteristics simulation results are post-processed and compared with the TLM calculation results. Fig.10 The experimental platform shown measures the dynamic characteristics of aerospace electromagnetic relays. First, remove the aerospace electromagnetic relay from its shell and fix it on the measuring table with a clamp, and aim the high-speed camera at the armature area of ​​the aerospace electromagnetic relay; then, open the current clamp and pass the coil wire through the opening of the current clamp, making sure that the wire is in the induction loop of the current clamp; finally, connect the high-speed camera to the computer, and the current clamp to the oscilloscope. At the moment when the aerospace electromagnetic relay is powered on, the high-speed camera will record the rotation process of the armature and display it on the computer, while the real-time changes of the coil current will be displayed and recorded on the oscilloscope. The dynamic characteristics of the experimental measurements are plotted on Fig.11 middle.

[0255] from Fig.11 It can be seen that the TLM numerical results of the present invention are consistent with COMSOL The dynamic characteristics results of the experimental measurements are very consistent, which proves the effectiveness of the three-dimensional digital prototype model of the electromagnetic-thermal-mechanical multi-physics field coupling of the aerospace electromagnetic relay proposed in the present invention. In the dynamic characteristics of the aerospace electromagnetic relay, when the switch of the control circuit is actuated, the coil begins to pass current and generate a magnetic field. In the initial stage, the mechanical resistance torque is greater than the electromagnetic torque, and the armature remains stationary. As the current gradually increases, the electromagnetic torque exerted on the armature will be greater than the mechanical resistance torque, and the armature will rotate. After a certain period of time, the armature will stop rotating after being attracted, the coil magnetic flux will rise to the saturation value, and the coil current will also tend to a steady state. These characteristics reveal the complex behavior and response characteristics of the aerospace electromagnetic relay during dynamic operation.

[0256] (2) Effect of temperature on dynamic characteristics

[0257] The temperature rise of aerospace electromagnetic relays during repeated actions will cause nonlinear changes in material properties, thereby affecting the dynamic characteristics of aerospace electromagnetic relays. The temperature rise experiment of aerospace electromagnetic relays is carried out under ambient room temperature. The temperature rise experiment test platform is as follows: Fig.12 As shown. Place the aerospace electromagnetic relay with the shell on the test bench and fix it with a fixture. Use the power-on mode of 1 second on and 1 second off (this mode is one of the typical working conditions for the electrical life test of aerospace electrical products) to continuously power on for 20 minutes. Use a temperature recorder to measure the temperature changes of the coil and the shell, and use a waveform recorder to record the changes of the coil current of the aerospace electromagnetic relay over time. Fig.13 The TLM and COMSOL simulations of aerospace electromagnetic relays working continuously for 20 minutes under the condition of power-on for 1 second and power-off for 1 second are given. and experimental measurement results.

[0258] observe Fig.13 It can be seen that the maximum temperature rise of the coil and shell of the aerospace electromagnetic relay after repeated operation for 20 minutes is 40K and 33K respectively. The coil current shows an overall downward trend as the internal temperature of the relay increases, and the coil current drops from 0.048A to 0.041A. This is because the temperature increase causes the resistance to increase. Fig.14 The temperature distribution shown in the figure shows that the internal temperature of the aerospace electromagnetic relay rises significantly due to the sealed shell, and the temperature rise is mainly concentrated in the coil and core area. In addition, in order to more intuitively compare the impact of temperature rise on dynamic characteristics, the shell of the aerospace electromagnetic relay was immediately removed after repeated action for 20 minutes, and the changes in dynamic characteristics were recorded with a high-speed camera and an oscilloscope. Fig.15 The dynamic characteristic curve of aerospace electromagnetic relay after repeated action for 20 minutes is given and compared with the dynamic characteristic before the action starts to analyze the influence of temperature rise on its dynamic characteristics.

[0259] from Fig.15 It can be seen from the dynamic characteristics displayed that the temperature rise after 20 minutes of repeated action will affect the dynamic characteristics of the aerospace electromagnetic relay. The coil current and electromagnetic torque become smaller due to the temperature rise, and the armature's pull-in time is also significantly prolonged. Therefore, during the reliability test and electrical life test of the aerospace electromagnetic relay, since the aerospace electromagnetic relay is in a state of long-term power-on and repeated opening and closing, the influence of temperature on the performance parameters of the electrical appliance cannot be ignored, which may affect its reliability and service life. Compared with the finite element model of a single physical field, the three-dimensional finite element model of multi-physical field coupling is more reliable and reasonable in practical engineering such as the design and reliability test of electrical products. In addition, the relative error between the calculation results of TLM and the experimental measurement value is controlled at about 4%, which verifies the accuracy and effectiveness of the three-dimensional digital prototype model proposed in the present invention under the repeated action conditions of electrical appliances.

[0260] (3) Comparison of solution efficiency

[0261] Since the number of mesh nodes and units of the three-dimensional FEM model of aerospace electromagnetic relay is large, in order to accelerate the convergence speed of the three-dimensional nonlinear problem of multi-physical field coupling, the present invention combines TLM and COMSOL The convergence error of is set to 1e-3. TLM and COMSOL Solve the dynamic characteristics of aerospace electromagnetic relays under the same computing environment. Fig.16 The dynamic characteristics solution information of two different scale grids solved by TLM with different numbers of threads is given in detail.

[0262] from Fig.16 It can be seen that the single-step solution time of TLM of the aerospace electromagnetic relay increases significantly when the armature is close to the energized position. This is because at the moment when the aerospace electromagnetic relay is energized, parameters such as the electromagnetic field and current density inside the relay will change suddenly, and the nonlinear effect will become more significant. At the same time, there is a small air gap at the energized position of the armature, which leads to an increase in the number of grids at the air gap and poor convergence of TLM. Therefore, the calculation efficiency of TLM at the armature energized position is greatly reduced. When the armature is energized, the coil current gradually tends to a steady state, the magnetic field no longer changes significantly, and the calculation time of TLM gradually decreases. In addition, from Fig.16 It can also be found that with the increase in the number of threads, the advantages of parallel computing can be brought into play, and the calculation time of TLM at each time step is significantly reduced. However, under the condition of coarse grid, the difference in solution time between 24 threads and 32 threads is not large, which means that when solving the digital prototype model of three-dimensional coarse grid, the acceleration effect achieved by continuing to increase the thread scale is not obvious. This is because with the increase in the number of threads, the allocation and coordination of parallel computing tasks may become more complicated, and too many threads will cause the communication and synchronization between threads to consume more time, which may offset the performance improvement brought by parallel computing. When the scale and complexity of the problem to be solved are further improved, such as solving the three-dimensional digital prototype model under fine grid conditions in the present invention, continuing to increase the number of threads can make full use of the parallel computing power provided by more threads to further improve the solution efficiency. Therefore, when using TLM to solve large-scale nonlinear problems, improving the parallel computing power of hardware can realize the accelerated calculation of digital prototype models, which is very suitable for large-scale network calculations.

[0263] In order to compare the parallel acceleration effect of TLM, Table 1 shows the TLM and COMSOL The total calculation time for solving the dynamic characteristics of aerospace electromagnetic relays under different meshes and different numbers of threads. As can be seen from Table 1, the increase in the number of meshes will lead to the increase in the number of threads in TLM and COMSOL. The solution time is significantly increased, but the calculation speed of TLM under multi-threaded parallel computing conditions is significantly faster than COMSOL Under coarse grid conditions, TLM parallel computing speed is at most faster than COMSOL It is 2.14 times faster, and the acceleration effect is not significant with the increase of the number of threads. This is because TLM takes a certain amount of time to allocate memory and interact with parallel data. For three-dimensional problems with a small number of grids, the increase in the number of threads cannot continue to significantly improve the calculation speed, and the calculation efficiency reaches saturation. As the number of grids increases, taking fine grids as an example, the parallel advantage of TLM gradually emerges, and the calculation time of TLM is much shorter than that of COMSOL. When using 16 threads, the computation time of TLM is 11.26 hours, and COMSOL The calculation time is 28.49 hours. The calculation speed of TLM is faster than that of COMSOL. 2.53 times faster. When the number of threads is 24, the TLM calculation time is 9.56 hours, COMSOL The calculation time is 25.33 hours, and the speedup ratio is 2.65, which is a significant acceleration effect. As the number of threads continues to increase, the speedup ratio reaches a maximum of 2.71 at 32 threads. Although the speedup ratio reaches the maximum at this time, the improvement rate of the speedup ratio is no longer significant, only 0.06 more than 24 threads. At this time, continuing to increase the number of threads will not achieve better parallel effects. Therefore, when using 24 threads, the parallel efficiency of TLM is higher and the parallel acceleration effect is better.

[0264] Table 1 Comparison of total calculation time and number of threads for three-dimensional dynamic characteristics of aerospace electromagnetic relays

[0265]

[0266] In summary, compared with the commercial finite element software COMSOL In comparison, under the condition of the same solution accuracy, the calculation speed of the three-dimensional TLM parallel algorithm in the case of the present invention is increased by 2.71 times at most. The improvement of the solution efficiency of the three-dimensional digital prototype model is of great significance for reducing the production cycle and production cost of aerospace electromagnetic relays. At the same time, TLM is very suitable for solving three-dimensional large-scale nonlinear problems, and its acceleration effect will increase with the improvement of the parallel capability of computer hardware, and has broad application prospects.

[0267] 3. Calculate the benefits

[0268] Compared with the single-field digital prototype model, the three-dimensional multi-physics field coupled digital prototype model based on TLM can comprehensively consider the influence of temperature on the dynamic characteristics of aerospace electromagnetic relays, and is compatible with the commercial finite element software COMSOL The relative error with experimental measurement is controlled within 4%. The established three-dimensional digital prototype model can accurately and quickly evaluate the dynamic characteristics and internal temperature changes of aerospace electromagnetic relays, proving the correctness and effectiveness of the three-dimensional multi-physics field coupled finite element parallel acceleration algorithm for aerospace electromagnetic relays proposed in the present invention.

[0269] Compared with the commercial finite element software COMSOL The TLM parallel algorithm has a good acceleration effect in solving three-dimensional multi-physics coupling problems, and can increase the calculation speed by up to 2.71 times. With the improvement of computer hardware parallel computing capabilities, the efficiency of TLM in solving three-dimensional multi-physics coupling problems can be further improved, and it has good application prospects.

[0270] The research results of the present invention have important theoretical significance for improving the environmental resistance and reliability design of aerospace electromagnetic relays, and also have reference significance for shortening the structural design and parameter optimization cycle of aerospace electromagnetic relays.

Claims

1. A three-dimensional electromagnetic relay multi-field coupling parallel calculation method based on the transmission line method, characterized in that The method comprises the following steps: Step 1: Establish a three-dimensional electromagnetic-thermal-mechanical multi-field coupling digital prototype model of aerospace electromagnetic relay: Step 1: Establish a three-dimensional finite element model of the electromagnetic system of the aerospace electromagnetic relay according to its geometric dimensions, discretize the solution domain into a finite number of tetrahedral units connected in size and shape, and control the number of grids according to the need for solution accuracy to complete the meshing of the finite element model; Step 1 and 2: Export and analyze the mesh information of the three-dimensional finite element model to obtain the number of tetrahedral elements, the number of nodes and the area to which the elements belong, as well as the node coordinates and node number mesh information of the tetrahedral elements, and calculate the edge lengths of all tetrahedral elements; Step 13: When the coil of the aerospace electromagnetic relay is powered on, a magnetic field will be generated. The transient magnetic field is described by the vector magnetic potential A, and its control equation is: Where A is the vector magnetic potential; ν is the magnetic resistivity; σ is the electrical conductivity; J e is the coil current density; t is the time; Step 14: Use edge finite element method to solve the three-dimensional finite element model: Step 141: Split the entire three-dimensional solution region Ω into many smaller tetrahedral units and the edges between tetrahedrons. e The vector magnetic potential vector A of any point in e (x,y,z), approximated by a linear combination of known interpolation functions: Where i is the number of edges of the tetrahedral unit; W i (x, y, z) is the interpolation function based on the edge vector; is the line integral of the vector magnetic potential field A(x, y, z) along the direction of the i-th edge; the interpolation function W i It is determined only by the four node position coordinates of the tetrahedral element, W i The projection on edge i is 1, and the projection on all other edges is 0; Step 142: The edge vector shape function W corresponding to the kth (k = 1, 2, ..., 6)th edge of the tetrahedral element k And its curl is expressed as: In the formula, l ij represents the edge length between node i and node j, and its positive and negative is determined by the direction vector; N1, N2, N3 and N4 represent the volume coordinates of the tetrahedral unit; and represents the gradient of the volume coordinates of the tetrahedral unit; Step 143. The gradient of the volume coordinates of the tetrahedral unit is expressed as: Where V e is the volume of the tetrahedral unit; and is the edge vector of the six edges of the tetrahedral unit; Step 144: For a single tetrahedral unit Ω e , whose four nodes are numbered n = 1, 2, 3, 4, then the volume V of the tetrahedral unit is e Calculated by the following formula: In the formula, x n ,y n ,z n They are the x, y, and z coordinates corresponding to the four nodes in the tetrahedral unit; Step 1, 4, and 5: Use the Galerkin method to interpolate the residual and the six interpolation functions W i Integrate the product and force the integral to be zero. After considering the Coulomb specification, the unit discrete equation in the form of a 10×10 matrix is ​​obtained: in Where N i is the volume coordinate of the tetrahedral unit; J e is the coil current density, the coil current is obtained through field-circuit coupling; A j is the line integral of the vector magnetic potential field A along the direction of the jth edge; Step 146: Establish a three-dimensional transient magnetic field transmission line model to solve the magnetic field equation in step 145. The specific steps are as follows: The transmission line method decouples the three-dimensional nonlinear system into a linear network and a nonlinear network, and obtains the equivalent circuit of the transmission line model. The equivalent circuit contains the circuit models of the incident stage and the reflection stage. The two stages are performed alternately in the linear network and the nonlinear network until convergence, thereby completing the solution of the three-dimensional nonlinear system, where: (1) The submatrix [K ij ] 6×6 and [D ij ] 6×6 The nonlinear resistance and linear capacitance of the equivalent circuit model are expressed as: (2) The nonlinear resistor is replaced by the transmission line and the reflection box in the reflection stage. The characteristic conductance of the transmission line is expressed by the following formula: (3) The linear capacitor is connected in parallel with the nonlinear resistor. The time derivative of the electromagnetic field is related to the capacitor. The capacitor is replaced by the resistor and the current source at the previous moment through the backward Euler method. The calculation formulas for admittance and current source are: Y Cij =C ij / Δt,Y Ci =C i / Δt,i,j∈[1,6] (4) In the incident phase, the node voltage A i and branch voltage (A i -A j ) enters the reflection box. In each reflection box, the incident voltage of one port will be scattered to all ports according to the coupling relationship of the nonlinear resistor. After the reflection stage, the incident voltage leaves the reflection box and becomes the reflected voltage. and reflected voltage The nonlinear relationship between is expressed as: The NR iteration is used to solve this 6×6 system of equations and the Jacobian matrix is ​​calculated according to the chain rule: After solving the 6×6 equations by NR iteration, the unit reluctance ν is obtained e The value of (5) In the incident phase, the reflected voltage containing the unit nonlinear resistance information is returned to the global circuit network, the transmission line and the reflected voltage are replaced by the Norton circuit, and then the node voltage of the global circuit network is updated through the linear network, and the incident voltage of the next iteration is calculated based on the node voltage. In addition, the circuit network in the incident phase is a large-scale linear network assembled by all unit admittances and Coulomb norms in the grid, and remains unchanged during the iteration process of each time step; Through steps one, four and six, the vector magnetic potential distribution of the magnetic field can be obtained in parallel; Step 147. According to the edge element formula, for each tetrahedral unit, the magnetic induction intensity B is represented by the vector magnetic potential A, and the three-dimensional nonlinear transient magnetic field equation is expressed as: Step 15: Use the field-circuit indirect coupling scheme to achieve the coupling between the relay electromagnetic field and the circuit: Step 151. The equation of the aerospace electromagnetic relay coil is expressed as: Where u is the coil voltage, i coil is the coil current, R is the coil resistance, ψ is the magnetic flux of the excitation coil; Step 152: Rewrite the formula in step 151 into an iterative form of a nonlinear equation: f (t) =i coil R+[ψ(t)-ψ(t-1)] / Δt-u=0 Where ψ(t) is the magnetic flux of the current time step, ψ(t-1) is the magnetic flux of the previous time step, and Δt is the discrete time step; Step 153: After obtaining the vector magnetic potential distribution of each edge of the coil unit according to the three-dimensional FEM model, the magnetic flux at the current time step is calculated by the following formula: Among them, N coil is the number of coil turns, S coil is the cross-sectional area of ​​the coil, n is the unit vector along the direction of the coil current, V coil is the volume of the coil area unit, N is the number of edges; Step 154: Use Newton's iteration method to construct the formula in step 152 into the following iterative format to calculate the current value i at the (t+1) moment: (t+1) , the iteration format is as follows: Get the current value i at the next moment (t+1) By inputting the transient magnetic field for iterative calculation, the vector magnetic potential value at the next moment is obtained, and the field-circuit coupling between the electromagnetic field of the aerospace electromagnetic relay and the control circuit can be realized; Step 16: Use the three-dimensional edge element and the unit load method to derive the parallel calculation method of electromagnetic force and electromagnetic torque: Step 161. According to the unit load method theory, the electromagnetic force is expressed as: Where q is the displacement direction of the armature, F q W is the electromagnetic force on the armature in the q direction. mag is the magnetic field energy; Step 162: In the tetrahedral unit Ω e The electromagnetic force on the midpoint k in the x direction is described as: Step 163: Combined with the three-dimensional nonlinear transient magnetic field equation, the unit node force F ek Simplified to: In the formula, Obtained through the BH curve of ferromagnetic materials, and Obtained by the following formula: Step 164. Volume V of tetrahedral unit e x k The partial derivative of is expressed as: Step 165: After calculating the electromagnetic force in the x direction, calculate the electromagnetic forces in the y and z directions in the same way; the electromagnetic force and electromagnetic torque are calculated in parallel; Step 17: Solve the electromagnetic field-mechanical motion coupling problem: Step 171. After obtaining the electromagnetic force of each node, calculate the electromagnetic torque of each node according to the following formula: Where M e is the electromagnetic torque of the e-th unit, F i e is the electromagnetic force magnitude of the node, is the vertical distance from the action point to the rotation axis, N n is the number of nodes in each unit; Step 172: Accumulate the electromagnetic torque of all contour nodes to obtain the total electromagnetic torque: Where M mag is the electromagnetic torque acting on the armature, N e is the number of tetrahedral units; Step 173. After calculating the electromagnetic torque, the dynamic equation of the armature of the aerospace electromagnetic relay in the rotational motion is described as: Where I is the moment of inertia, α is the angular acceleration of the armature, θ is the rotation angle of the armature, and M r is the mechanical resistance torque acting on the armature; Step 18: Finite element calculation of thermal conduction of aerospace electromagnetic relay: Step 181. The three-dimensional partial differential equation of the transient temperature field is expressed as: In the formula, ρ is the density of the material, C is the specific heat capacity of the material, T is the temperature of the solid or air, λ is the thermal conductivity of the medium, and Q is the heat generation power per unit volume; Step 182: Use the tetrahedron node temperature distribution to characterize the temperature field, and use the FEM discretization process to obtain the following equation: After the above formula is discretized in time using the backward Euler method, the following algebraic equation is obtained: In the formula, Where l,n=1,2,3,4; is the density of the material in the unit area, C e is the specific heat capacity of the material in the unit area, λ e is the thermal conductivity of the unit medium, G e is the internal heat source intensity, p l ,q l 、r l and l is a function related to the node coordinates x, y and z; Step 183. The heat sources in aerospace electromagnetic relays mainly include Joule loss and eddy current loss, of which eddy current loss P eddy After the magnetic field calculation is completed, it is obtained through post-processing and expressed as: Where E is the electric field strength, Ω core It is the ferromagnetic region; Step 2: Solve the electromagnetic-thermal-mechanical multi-field coupled digital prototype model of the three-dimensional finite element model of the aerospace electromagnetic relay: Step 21: Read the temperature distribution at the last moment and update the material properties of each area; Step 22: Read the vector magnetic potential distribution and magnetic reluctance at the current moment, and calculate the Joule loss and eddy current loss; Step 2: Calculate the temperature field to update the temperature distribution, and update the conductivity of the coil and the ferromagnetic material; Step 24: Calculate the electromagnetic force and electromagnetic torque on the armature according to the unit load method; Step 25: Calculate the mechanical motion equation to determine whether the armature rotates. If so, update the mesh information of the three-dimensional electromagnetic-thermal-mechanical multi-field coupled digital prototype model. If the armature remains stationary, there is no need to update the mesh information. Step 26: Calculate the coil current value of the next time step based on the field-circuit coupling, and save the vector magnetic potential and magnetic reluctance of the next time step.

2. The three-dimensional electromagnetic relay multi-field coupling parallel calculation method based on the transmission line method according to claim 1 is characterized in that The W mag The calculation formula is expressed as: Where W e is the unit magnetic field energy, N e is the number of tetrahedral elements that include nodal forces.

3. The three-dimensional electromagnetic relay multi-field coupling parallel calculation method based on the transmission line method according to claim 1 is characterized in that The p l ,q l 、r l and l The functional relationship is as follows:

Images