A precise levitation force mechanism modeling method for spherical magnetic levitation flywheel motor

Abstract

The present invention provides a method for accurately modeling the levitation force mechanism of a spherical magnetic levitation flywheel motor. The method comprises the following steps: Step 1: obtaining the radial and axial equivalent magnetic circuits of the levitation poles of the spherical magnetic levitation flywheel motor based on the equivalent magnetic circuit method; Step 2: analyzing the diffuse magnetic flux of the levitation poles of the spherical magnetic levitation flywheel motor using finite element analysis; Step 3: using the equivalent magnetic circuit to obtain the air gap permeance that accounts for rotor eccentricity, and then correcting the air gap main flux density to obtain the corrected air gap main flux density; Step 4: analyzing different rotor position angle ranges based on the equivalent magnetic circuit model of the air gap permeance, combining the air gap length and the air gap main flux density, calculating the air gap permeance of each part, and establishing a levitation force model for the spherical magnetic levitation flywheel motor. This modeling method reduces the modeling error from 20% to 5%, providing an effective method and model for accurately characterizing the dynamic characteristics of the motor under rotor eccentricity, magnetic circuit saturation, and diffuse leakage, laying a foundation for the analysis, design, and operation control of spherical magnetic levitation flywheel motors.

Description

Technical Field

[0001] The present invention relates to the technical field of magnetic levitation motors, and in particular to a method for modeling the precise levitation force mechanism of a spherical magnetic levitation flywheel motor. Background Art

[0002] Flywheel energy storage systems (FESS) have attracted considerable attention due to their high power density, long cycle life, fast charge and discharge, and clean, environmentally friendly design. They hold broad application prospects in uninterruptible power supplies, electric vehicles, aerospace, and other fields. Magnetic levitation motors, combining the dual functions of magnetic bearing suspension and motor rotation, offer advantages such as wear-free and lubrication-free operation, high power, and ultra-high-speed operation. Incorporating a magnetic levitation flywheel machine (BFM) into a FESS reduces system size and losses, while increasing critical speed and power density. This makes it an ideal choice for FESS suspension support and energy conversion.

[0003] The BFM mechanism model is a prerequisite for its optimized design, performance analysis, and efficiency evaluation, and also lays the foundation for high-performance control, attracting significant attention from scholars both domestically and internationally. In the 1990s, Japanese scholar Chiba A et al. developed a mathematical model for the 12 / 8 dual-winding BFM structure. In China, Deng Zhiquan et al. at Nanjing University of Aeronautics and Astronautics proposed a mathematical model based on the air gap permeance obtained from a linear magnetic circuit and a modified elliptical magnetic circuit. Sun Yukun et al. at Jiangsu University further expanded on this approach. Ge Baoming et al. at Beijing University of Aeronautics and Astronautics developed a suspension force and torque model based on an analysis of its electromagnetic properties.

[0004] Existing literature generally uses the principle of virtual work to derive models. By obtaining the inductance matrix, the virtual displacement method is used to derive the suspension force model. This method accounts for the effects of eccentricity to some extent but neglects the effects of magnetic circuit saturation. Finite element analysis can account for magnetic saturation, but its computational speed is slow, making it difficult to apply in practical control systems. The Maxwell tensor method, which takes into account magnetic saturation characteristics, is fast and concise, but requires manual empirical correction. Model accuracy, particularly under conditions of rotor eccentricity, magnetic circuit saturation, and diffuse flux, needs to be improved. Summary of the Invention

[0005] In response to the model mismatch problem of existing mechanism models under rotor eccentricity, magnetic circuit saturation and diffuse magnetic flux, the present invention provides a method for modeling the precise suspension force mechanism of a spherical magnetic levitation flywheel motor. Based on the analysis of the working principle and magnetic flux distribution characteristics of the spherical magnetic levitation flywheel motor, the present invention adopts a magnetic field segmentation method to solve the magnetic permeability of different diffuse magnetic flux regions of the motor. According to the inherent characteristics of the flywheel motor rotor eccentricity and magnetic circuit saturation, the Maxwell tensor method and the equivalent magnetic circuit method are combined to construct a precise suspension force model of the motor. The results of finite element simulation and comparative analysis show that the error of the proposed model is reduced from 20% to 5%, which verifies that the constructed model more accurately describes the dynamic performance under rotor eccentricity, magnetic circuit saturation and diffuse magnetic flux.

[0006] In order to achieve the above purpose, the present invention is achieved through the following technical solutions:

[0007] The present invention is a method for accurately modeling the levitation force mechanism of a spherical magnetic levitation flywheel motor, comprising the following steps:

[0008] S1. Obtain the equivalent magnetic circuits of the radial and axial directions of the spherical magnetic levitation flywheel motor's suspension poles based on the equivalent magnetic circuit method.

[0009] S2. Use finite element analysis to analyze the diffusion flux of the suspended pole of the spherical magnetic levitation flywheel motor. In the present invention, within a 30° cycle of rotor movement, based on the spherical magnetic field distribution law and the magnetic flux path, the rotor position angle range is divided into two parts, and the magnetic field segmentation method is used to divide the entire air gap magnetic field into several flux pipes. The magnetic permeability of each flux pipe is calculated separately, and then the magnetic flux pipes are synthesized to calculate the entire air gap magnetic permeability. Analysis in the range of rotor position angle -7.5°<θ<7.5° shows that due to the symmetry of the magnetic poles, the diffuse magnetic flux is symmetrical about the two sides of the magnetic poles. The magnetic field is divided on one side and divided into region 1, the extreme end face part of the magnetic pole, that is, the flux pipe of the main magnetic flux between the poles, and the corresponding magnetic permeability is G ax1 Region 2 is a 1 / 4 solid cylinder, which is the flux conduit for the diffusion flux, and the corresponding magnetic permeance is G s1 ; Region 3 is a 1 / 4 hollow cylinder, which is the flux conduit for the diffusion flux, and the corresponding magnetic permeance is G s2 The analysis of the diffusion flux of other suspended poles is similar, and the equivalent magnetic circuit of the spherical magnetic levitation flywheel motor considering the diffusion flux is obtained;

[0010] S3. Establishing an equivalent magnetic circuit of a spherical magnetic levitation flywheel motor taking into account the diffuse magnetic flux using step S2. Analyzing and calculating the change in the air gap length between the rotor teeth and the suspension pole shoe teeth due to the rotor eccentricity, obtaining the air gap permeability taking into account the rotor eccentricity. Furthermore, considering the need for magnetic circuit saturation, correcting the air gap main flux density to obtain the corrected air gap main flux density.

[0011] S4. Based on the equivalent magnetic circuit model of the air gap permeability obtained in step S2, combined with the air gap length and air gap main flux density obtained in step S3, analyze the conditions of different rotor position angle ranges, calculate the air gap permeability of magnetic flux ducts of different shapes after considering rotor eccentricity, and establish a spherical magnetic levitation flywheel motor suspension force model.

[0012] Preferably, in step S1, the magnetic field distribution of the spherical magnetic levitation flywheel motor is obtained by finite element analysis, and the radial and axial equivalent magnetic circuits are obtained based on the equivalent magnetic circuit method.

[0013] Preferably, in step S2, the rotor rotates for a period of 30°.

[0014] Preferably, in step S2, one rotation cycle of the rotor is divided into two parts: -7.5°<θ<7.5°, -15°<θ<-7.5°, and 7.5°<θ<15°.

[0015] The segmented magnetic field method is to divide the entire air gap flux, including the edge flux, into several flux tubes with simple geometric shapes according to their possible paths. The permeance of each flux tube is first calculated separately, and then the permeances of the parallel flux tubes are added together to obtain the total air gap flux.

[0016] Preferably, in step S2, the air gap magnetic circuit is divided into magnetic flux tubes of different shapes. The shapes of the magnetic flux tubes are rectangular magnetic flux tubes under the magnetic poles, 1 / 4 solid cylinders, and 1 / 4 hollow cylinders. The magnetic permeance calculation formulas of the three shapes of magnetic flux tubes are:

[0017] The permeance of the rectangular flux tube under the magnetic pole is:

[0018]

[0019] Magnetic permeance of a 1 / 4 solid cylindrical flux tube:

[0020]

[0021] Magnetic permeance of a 1 / 4 hollow cylindrical flux tube:

[0022]

[0023] Where μ0 is the vacuum magnetic permeability, η is the arc length of the rotor tooth corresponding to the sphere, R is the radius of the rotor tooth corresponding to the sphere, h1 is the axial height of the spherical rotor tooth, δ av1 is the average air gap length of region 1, is the angle between the axial rotor tooth spherical surface, A1 is the average cross-sectional area of ​​the rectangular flux tube in region 1, V is the volume of 1 / 4 solid cylindrical flux tube, A2 is the average cross-sectional area of ​​1 / 4 hollow cylindrical flux tube, δ av2 is the average magnetic flux length in region 2, δ av3 is the average magnetic flux length in region 3.

[0024] Preferably, in step S3, the influence of rotor eccentricity on the air gap magnetic flux is considered, and the relationship of the air gap length after the rotor is eccentric bmm in the negative direction of the y-axis is:

[0025]

[0026] Further deduction:

[0027] δ1≈δ0+bsinθ m

[0028] When the rotor is offset a mm in the positive direction of the x-axis, the air gap length relationship is:

[0029]

[0030] Further deduction:

[0031] δ1≈δ0+acosθ m

[0032] When the rotor is offset a mm in the positive direction of the x-axis and b mm in the negative direction of the y-axis, the relationship between the air gap lengths is:

[0033]

[0034] Further deduction:

[0035] δ1≈δ0+acosθ m -bsinθ m

[0036] Among them, δ0 is the average air gap length when the rotor is not eccentric, δ1 is the air gap length between the stator and rotor poles after eccentricity, and θ m is the mechanical position angle, and r is the radius of the sphere where the stator is located.

[0037] Preferably, the main magnetic flux density of the air gap in step S3 is:

[0038]

[0039] The air gap diffusion flux density is:

[0040]

[0041] Combined with the finite element analysis results, the air gap magnetic flux correction formula is designed as follows:

[0042]

[0043] Where N is the number of turns of the suspension winding, I0 is the current provided by the permanent magnet, I c is the control current of the suspension winding; B nlis the air gap flux density after correction; B c B is the critical value of the magnetic density of ferromagnetic materials tending to saturation; l According to formula B m 、B q The air gap magnetic flux density is obtained; k is the correction coefficient.

[0044] Preferably, in step S4, when the rotor position angle is -7.5°<θ<7.5°, the radial spherical air gap reluctance is:

[0045]

[0046] The magnetic flux generated by the ±x-direction branch and the ±y-direction branch is:

[0047]

[0048] Among them, G zax1 , G zax2 , G zay1 , G zay2 They are the total magnetic permeability of the air gap in the ±x and ±y directions of phase A, Other similarities; N is the number of turns of the A-phase suspension winding, I0 is the current provided by the permanent magnet, I c is the control current of the A-phase suspension winding, C r is the magnetic flux leakage coefficient of the radial part of the motor, Φ x1 , Φ x2 They are the magnetic flux generated by the ±x direction branches, Φ y1 , Φ y2 They are the magnetic flux generated by the ±x direction branches, and the unit is Weber (Wb).

[0049] Preferably, the electromagnetic force F exerted on the rotor in step S4 is x for:

[0050]

[0051] Wherein, α is half of the spherical angle of the suspended pole, and A represents the area of ​​the pole end surface including the diffusion flux range.

[0052] Preferably, in step S4, when the rotor position angle is 15°<θ<-7.5° and 7.5°<θ<15°, the irregular diffusion magnetic flux at both ends is taken as a semicircular integral path dL and the electromagnetic force dF is decomposed to obtain the radial force dF r and tangential force dF t for:

[0053]

[0054] The radial and tangential components of the electromagnetic force generated on this path are decomposed into:

[0055]

[0056] When n approaches infinity, the radial and tangential components are:

[0057]

[0058] in, h is the thickness of the suspension of phase A or phase B, i is the number of parts, l 78 、l 56 is the integral path length; μ0 is the vacuum permeability.

[0059] Preferably, the levitation force generated by the irregular diffusion magnetic flux in step S4 is:

[0060]

[0061]

[0062] in, B f , B f1 , B x1f1 , B x2f1 , B x1f , B x2f According to B q Formula for the edge air gap flux density after correction in ±x direction.

[0063] The beneficial effects of the present invention are:

[0064] 1. Based on the principle of minimum magnetic resistance, the present invention divides a rotor rotation cycle of 30° into two parts: -7.5°<θ<7.5°, -15°<θ<-7.5°, and 7.5°<θ<15°, based on the diffusion flux path. This achieves accurate derivation of the air gap permeability and improves the accuracy of the model.

[0065] 2. The present invention takes into account both rotor eccentricity and diffuse magnetic flux, and combines the air gap flux correction formula to describe the magnetic circuit saturation characteristics. The finite element simulation comparative analysis of the obtained model shows that the error of the model built by the proposed modeling method is reduced from 20% to 5%, which provides an effective method, model and idea for accurately characterizing the dynamic characteristics of the motor under rotor eccentricity, magnetic circuit saturation and diffuse leakage, and lays a theoretical foundation for the analysis, design and operation control of spherical magnetic levitation flywheel motors. BRIEF DESCRIPTION OF THE DRAWINGS

[0066] Figure 1 This is a flow chart of the precise suspension force mechanism modeling method for a spherical magnetic levitation flywheel motor proposed in the present invention.

[0067] Figure 2This is the single-phase magnetic pole magnetic field distribution diagram of the spherical magnetic levitation flywheel motor proposed in the present invention.

[0068] Figure 3 This is a schematic diagram of the magnetic flux pipeline division of the spherical magnetic levitation flywheel motor proposed in the present invention.

[0069] Figure 4 It is a schematic diagram of the air gap integral path of the spherical magnetic levitation flywheel motor proposed in the present invention.

[0070] Figure 5 It is the equivalent magnetic circuit of the spherical magnetic levitation flywheel motor considering the diffusion magnetic flux.

[0071] Figure 6 It is a schematic diagram of the air gap length model of the rotor offset in the negative direction of the y-axis.

[0072] Figure 7 It is a schematic diagram of the rotor air gap length model offset in the positive direction of the x-axis.

[0073] Figure 8 This is a schematic diagram of the rotor air gap length model that is simultaneously offset in the positive direction of the x-axis and the negative direction of the y-axis.

[0074] Figure 9 This is the radial and axial equivalent magnetic circuit diagram of the spherical magnetic levitation flywheel motor proposed in the present invention.

[0075] Figure 10 This is a schematic diagram of electromagnetic force calculation for air gap integral paths 7 to 8.

[0076] Figure 11 These are the structural parameters of the spherical magnetic levitation flywheel motor proposed in the present invention.

[0077] Figure 12 It is a schematic structural diagram of the three-dimensional finite element model of the spherical magnetic levitation flywheel motor proposed in the present invention.

[0078] Figure 13 This is a graph showing the relationship between the suspension force on the motor rotor, the x-direction offset distance, and the mechanical position angle.

[0079] Figure 14 This is a graph showing the relationship between the suspension force on the motor rotor, the y-direction offset distance, and the mechanical position angle.

[0080] Figure 15 This is a graph showing the relationship between the levitation force of a spherical magnetic levitation flywheel motor and the ampere-turns of the levitation winding.

[0081] Figure 16 This is the relationship diagram between the motor's magnetic saturation model and the suspension force obtained by simulation and the winding ampere-turns. DETAILED DESCRIPTION

[0082] The following drawings illustrate embodiments of the present invention. For clarity, many practical details are included in the following description. However, it should be understood that these practical details are not intended to limit the present invention. In other words, in some embodiments of the present invention, these practical details are not essential. Furthermore, to simplify the drawings, some commonly used structures and components are depicted in simplified schematic form.

[0083] like Figure 1 As shown, the present invention is a method for modeling the precise levitation force mechanism of a spherical magnetic levitation flywheel motor, and the modeling method includes the following steps:

[0084] Step 1: Obtain the magnetic field distribution of the spherical magnetic levitation flywheel motor by finite element analysis, and obtain the equivalent magnetic circuits of the radial and axial directions of the levitation poles of the spherical magnetic levitation flywheel motor based on the equivalent magnetic circuit method.

[0085] Step 2: Use finite element analysis to analyze the diffusion flux of the suspended pole of the spherical magnetic levitation flywheel motor. In the present invention, within a 30° period of rotor movement, based on the spherical magnetic field distribution law and the magnetic flux path, the rotor position angle range is divided into two parts, and the magnetic field segmentation method is used to divide the entire air gap magnetic field into several flux pipes. The magnetic permeability of each flux pipe is calculated separately, and then the magnetic flux pipes are synthesized to calculate the entire air gap magnetic permeability. The rotor position angle is analyzed in the range of -7.5°<θ<7.5°. Taking a suspended pole on the x-axis of phase A as an example, due to the symmetry of the magnetic pole, the diffused magnetic flux is symmetrical about the two sides of the magnetic pole. The magnetic field is divided on one side and divided into region 1, the extreme end face part of the magnetic pole, that is, the flux pipe of the main magnetic flux between the poles. The corresponding magnetic permeability is G ax1 Region 2 is a 1 / 4 solid cylinder, which is the flux conduit for the diffusion flux, and the corresponding magnetic permeance is G s1 ; Region 3 is a 1 / 4 hollow cylinder, which is the flux conduit for the diffusion flux, and the corresponding magnetic permeance is G s2 The diffuse magnetic flux on the other side of the suspension pole and at both ends of the suspension pole on the x and y axes are also divided into three regions to calculate the corresponding magnetic permeability, thereby obtaining the equivalent magnetic circuit of the spherical magnetic levitation flywheel motor considering the diffuse magnetic flux.

[0086] Step 3: Use step 2 to establish an equivalent magnetic circuit of a spherical magnetic levitation flywheel motor taking into account the diffuse magnetic flux. By analyzing and calculating the change in the air gap length between the rotor teeth and the suspended pole shoe teeth due to rotor eccentricity, the air gap magnetic permeability considering the rotor eccentricity is obtained. In addition, considering the need for magnetic circuit saturation, the air gap main magnetic flux density is corrected to obtain the corrected air gap main magnetic flux density.

[0087] Step 4: Based on the equivalent magnetic circuit model of the air gap permeability obtained in step 2, combined with the air gap length and air gap main flux density obtained in step 3, analyze the conditions of different rotor position angle ranges, calculate the air gap permeability of each part, and establish the spherical magnetic levitation flywheel motor suspension force model.

[0088] In step 1, the structural parameters of the spherical motor are as shown in the attached Figure 11 As shown in the figure, the structural parameters of the spherical motor include the stator outer diameter, rotor outer diameter, rotor inner diameter, spherical rotor diameter, spherical stator diameter, permanent magnet inner diameter, permanent magnet outer diameter, permanent magnet thickness, air gap length, rotor pole arc, torque pole arc, suspension pole arc, suspension winding turns and axial length. The structural diagram of the 3D finite element model of the spherical magnetic levitation flywheel motor is shown in the attached figure. Figure 12 The magnetic field distribution of the spherical magnetic levitation flywheel motor is obtained by finite element analysis, and the radial and axial equivalent magnetic circuits are obtained based on the equivalent magnetic circuit method, as shown in the attached figure. Figure 9 As shown, F xa1 、F xa2 is the magnetomotive force of the A-phase suspension winding in the ±x direction, F ya1 、F ya2 is the magnetomotive force of the A-phase suspension winding in the ±y direction, F pm is the magnetomotive force of the permanent magnet; R xa1 、R xa2 is the magnetic reluctance of the suspended pole core of phase A in the ±x direction, R ya1 、R ya2 is the magnetic reluctance of the suspended pole core of phase A in the ±x direction, R xb1 、R xb2 is the magnetic reluctance of the B-phase suspended pole core in the ±x direction, R ry is the rotor core reluctance, R pm is the permanent magnet reluctance, R xaq1 、R xaq2 is the air gap reluctance of phase A in the ±x direction, R yaq1 、R yaq2 is the air gap reluctance of phase A in the ±y direction, R xbq1 、R xbq2 is the air gap reluctance of phase B in ±x direction.

[0089] In step 2, the rotor position angle θ is defined as the angle of the rotor rotating counterclockwise from the position where the suspension pole axis coincides with the rotor pole axis. That is, the initial position angle θ = 0° is the position where the suspension pole axis coincides with the rotor pole axis. The air gap magnetic field line distribution of the spherical magnetic levitation flywheel motor obtained by finite element analysis is shown in the attached figure. Figure 2 The figure shows two main components: the main magnetic flux between the poles and the diffuse magnetic flux at the edges. Based on the principle of minimum magnetic resistance, most of the magnetic flux flows from the suspended pole shoe tooth surface to the rotor tooth surface, while a small portion flows into the rotor tooth side surface and the rotor yoke. Based on the diffuse magnetic flux path, a 30° rotor rotation period is divided into two parts: -7.5° < θ < 7.5°; -15° < θ < -7.5°; and 7.5° < θ < 15°.

[0090] The analysis at -7.5°<θ<7.5° is as shown in the attached figure. Figure 3As shown in the three-dimensional structure diagram, due to the symmetry of the magnetic poles, the diffusion flux is symmetrical about the two sides of the magnetic poles. The magnetic field on one side is divided into three regions as shown in the figure. Region 1 is the end face of the magnetic pole, which is the flux conduit of the main magnetic flux between the poles. The corresponding magnetic permeance is G ax1 The calculation formula is shown in formula (1); Region 2 is a 1 / 4 solid cylinder, which is the flux conduit for the diffusion flux, and the corresponding magnetic permeance is G s1 The calculation formula is shown in formula (2); Region 3 is a 1 / 4 hollow cylinder, which is the flux conduit for the diffusion flux, and the corresponding magnetic permeance is G s2 The calculation formula is shown in formula (3).

[0091] When the rotor is at 15°<θ<-7.5° and 7.5°<θ<15°, the diffusion flux path will change and become asymmetric at both ends, as shown in the following figure. Figure 4 As shown in the figure, the calculation method of the air gap magnetic permeability of paths 1 to 5 is the same as that of -7.5°<θ<7.5°. The difference is that the air gap magnetic flux of paths 5 to 8, that is, the diffusion magnetic flux between the rotor pole and the suspended pole air gap, cannot be calculated according to the division of the magnetic flux pipeline, and the air gap magnetic flux path needs to be integrated to obtain it.

[0092]

[0093]

[0094]

[0095] Where μ0 is the vacuum magnetic permeability, η is the arc length of the rotor tooth corresponding to the sphere, R is the radius of the rotor tooth corresponding to the sphere, h1 is the axial height of the spherical rotor tooth, δ av1 is the average air gap length in region 1, β is the angle between the axial rotor tooth spherical surface, A1 is the average cross-sectional area of ​​the rectangular flux tube in region 1, V is the volume of 1 / 4 solid cylindrical flux tube, A2 is the average cross-sectional area of ​​1 / 4 hollow cylindrical flux tube, δ av2 is the average magnetic flux length in region 2, δ av3 is the average magnetic flux length in region 3.

[0096] In step 3, based on the actual structure of the spherical magnetic levitation flywheel motor, combined with step 1 to obtain the radial and axial equivalent magnetic circuits of phase A and phase B, and step 2 to consider the diffuse magnetic flux and use the magnetic field segmentation method to obtain the edge air gap magnetic permeability, the magnetic flux passes through the edge area of ​​the suspension pole and the end face of the magnetic pole at the same time, so the magnetic resistance of the two parts will be in parallel, and the magnetic permeability of the two parts will be in series. The permanent magnet whose magnetic permeability does not change during the movement of the rotor is equivalent to a flux tube, and the magnetic circuit is represented by magnetomotive force and magnetic permeability. After the control current is passed through the suspension pole, it is also represented by magnetomotive force. After the magnetic flux flows into the rotor from the air gap, the magnetic flux generated by the suspension winding forms a closed loop through the radial magnetic circuit, and the bias magnetic flux generated by the permanent magnet forms a closed loop through the axial magnetic circuit. As shown in the attached figure Figure 5 As shown in the figure, G ax1 ,G ax2 ,G ay1 ,G ay2 ,G bx1 ,G bx2 ,G by1 ,G by2 The air gap permeability under the extreme ends of the x-axis and y-axis magnetic poles of phase A and phase B, G s is the air gap permeability at the edge of the diffusion flux; Ni ax1 ,Ni ax2 ,Ni ay1 ,Ni ay2 ,Ni bx1 ,Ni bx2 ,Ni by1 ,Ni by2 is the magnetomotive force of the suspension winding; F pm is the magnetomotive force of the permanent magnet; G0 is the magnetic permeance of the permanent magnet.

[0097] Since the air gap is very small, its slight change will have a great impact on the radial force on the rotor, so it is necessary to consider the effect of rotor eccentricity on the air gap flux. Figure 6 As shown, the relationship between the air gap length after the rotor is eccentric b mm in the negative y-axis direction is shown in formula (4), and the result is shown in formula (5). The magnetic circuit in the rotor pole overlap area is set to a straight magnetic circuit.

[0098]

[0099] Further deduction:

[0100] δ1≈δ0+bsinθ m (5)

[0101] As attached Figure 7 The figure shows the schematic diagram of the air gap length model after the rotor is offset a mm in the positive direction of the x-axis. The relationship of the air gap length is shown in formula (6), and the result is shown in formula (7).

[0102]

[0103] Further deduction:

[0104] δ1≈δ0+acosθ m (7)

[0105] As attached Figure 8 The figure shows the schematic diagram of the air gap length model after the rotor is offset a mm in the positive direction of the x-axis and b mm in the negative direction of the y-axis. The relationship of the air gap length is shown in formula (8), and the result is shown in formula (9).

[0106]

[0107] Further deduction:

[0108] δ1≈δ0+acosθ m -bsinθ m (9)

[0109] Among them, δ0 is the average air gap length when the rotor is not eccentric, δ1 is the air gap length between the stator and rotor poles after eccentricity, and θ m is the mechanical position angle, and r is the radius of the sphere where the stator is located. From the above, it can be seen that δ1 is not only related to the rotor offset, but also to the mechanical position angle θ. m Because the overlap between the suspended pole and the rotor pole axis at the initial position of the rotor is θ m The range of change is -15°<θ m <15°, so sinθ m ≈θ m , cosθ m ≈1.

[0110] In actual situations, when the current increases to a certain level, the magnetic field inside the motor becomes saturated. At this time, the air gap flux density of the motor will no longer follow the linear calculation formula, as shown in the attached figure. Figure 4 The air gap main flux density (Equation (10)) and the air gap diffuse flux density (Equation (11)) are calculated and corrected. The finite element analysis uses the nonlinear ferromagnetic material DW360_50. Therefore, it is necessary to correct the formula for calculating the suspension force generated by the motor in the saturated state.

[0111]

[0112]

[0113] Where: N is the number of turns of the suspension winding, I0 is the current provided by the permanent magnet, I c is the control current of the suspension winding.

[0114] Combined with the finite element analysis results, the present invention designs the air gap magnetic flux correction formula:

[0115]

[0116] Where: N is the number of turns of the suspension winding, I0 is the current provided by the permanent magnet, I c is the control current of the suspension winding; B nl is the air gap flux density after correction; B c B is the critical value of the magnetic density of ferromagnetic materials tending to saturation; l According to formula B m 、B q The air gap magnetic flux density is obtained; k is the correction coefficient.

[0117] In step 4, using phase A of the spherical magnetic levitation flywheel motor as an example, the radial levitation force is calculated by combining the equivalent magnetic circuit of step 1, the air gap permeance of step 2, and the air gap length and air gap flux density of step 3. When the rotor position angle is -7.5° < θ < 7.5°, the radial spherical air gap reluctance is given by equation (13). The magnetic flux generated by the ±x- and ±y-direction branches is then calculated as shown in equation (14).

[0118]

[0119]

[0120] Where: G zax1 , G zax2 , G zay1 , G zay2 They are the total magnetic permeability of the air gap in the ±x and ±y directions of phase A, G zax2 , G zay1 , G zay2 Calculation method and G zax1 Same; N is the number of turns of the A-phase suspension winding, I0 is the current provided by the permanent magnet, I c is the control current of the A-phase suspension winding, C r is the magnetic flux leakage coefficient of the radial part of the motor, Φ x1 , Φ x2 They are the magnetic flux generated by the ±x direction branches, Φ y1 , Φ y2 are the magnetic flux generated by the ±x direction branches, and the unit is Weber (Wb). Combining formula (13) and formula (14), the electromagnetic force F on the rotor can be obtained. x As shown in formula (15).

[0121]

[0122] Where: α is half of the spherical angle of the suspended pole, and A represents the area of ​​the pole end surface including the diffusion flux range.

[0123] When the rotor position angle is 15°<θ<-7.5° and 7.5°<θ<15°, as shown in the attached Figure 4 The calculation method of the air gap magnetic permeability of the path 1 to 5 is the same as that of -7.5°<θ<7.5°, that is, the suspension force of this part is F x In addition, the suspension force generated by the irregular diffusion flux at the upper and lower ends must be derived. The integral path is a straight line, and the edge flux integral path 7 to 8 is approximately a semicircle (the symmetrical part is not involved in the calculation), as shown in the attached figure. Figure 10 As shown. Select the straight line where the radius is at 90° in the figure as the symmetry axis, take the integral length of the symmetric position dL of the symmetry axis on the arc, and decompose the electromagnetic force dF to obtain the electromagnetic force generated in the radial and tangential directions are equal, as shown in formula (16). Divide the arc where the 7-8 path is located into 2n parts, n parts on each side of the symmetry axis, and the length of each part is Therefore, the radial and tangential components of the electromagnetic force generated on the path 7 to 8 are decomposed as shown in formula (17). When n tends to infinity, formula (17) can be expressed as formula (18).

[0124]

[0125]

[0126]

[0127] Where, h is the thickness of the suspension of phase A or phase B, i is the number of parts, l 78 、l 56 is the integral path length; μ0 is the vacuum permeability.

[0128] Therefore, the suspension force generated by the irregular diffusion flux at the upper and lower ends is as shown in formula (19). Therefore, the total radial force generated by this part in the x-axis direction of phase A is as shown in formula (20).

[0129]

[0130]

[0131] Where: B f , B f1 , B x1f1 , B x2f1 , B x1f , B x2f The marginal air gap flux in the ±x direction after correction in step 3 is obtained according to formula (11) and formula (12) respectively. In summary, when the rotor position angle is -15°<θ<-7.5° and 7.5°<θ<15°, the suspension force expression is F=Fx +F r .

[0132] Simulation verification:

[0133] To further illustrate the effectiveness of the precise levitation force mechanism model of the spherical magnetic levitation flywheel motor constructed by the method of the present invention, the finite element analysis (FEA) of the motor and the model constructed by the present invention are compared in terms of controlling the current I c =0.5A, attached Figure 13 The relationship between the rotor's suspension force, the x-direction offset distance, and the mechanical position angle is shown in the figure below: Figure 14 The relationship between the rotor's suspension force, the y-direction offset distance, and the mechanical position angle is shown in Figure 2. The model result is slightly larger than the simulation value. This is because the eddy current effect is ignored before deriving the mathematical model. Since the error is not large and the outer rotor is not a solid iron core, the eddy current effect has little effect on the suspension force model and can be ignored, meeting the accuracy requirements of the modeling.

[0134] In order to further illustrate the effectiveness of the model, the present invention uses the established model to solve the relationship between the suspension force and the ampere-turns of the suspension winding as shown in the attached figure. Figure 15 As shown in the figure, compared with the FEA analysis results, the model results better reflect the relationship between the suspension force and the winding ampere-turns.

[0135] Attachment Figure 16 A graph shows the relationship between the linear model and the magnetic saturation model of the present invention, as well as the simulated levitation force as the winding ampere-turns vary from 0 to 600 AT. The levitation pole and torque pole are decoupled, and considering the rotor eccentricity, when current is applied to the winding in the x direction and no current is applied to the winding in the y direction, the levitation force curves are shown at rotor position angles θ = 0°, -2.5°, -5°, -7.5°, -10°, and -15°, respectively. It can be found that when the current is high and the motor experiences magnetic circuit saturation, the levitation force described by the linear model far exceeds the normal error limit and can no longer be used to describe the motor's levitation force characteristics. The calculated values ​​of the levitation force model proposed in the present invention are basically consistent with the FEA analysis results, and can more accurately describe the levitation force variation curve of the motor in the high current saturation state, which is consistent with the magnetic circuit saturation characteristics of the motor during actual operation.

[0136] This paper proposes a method for accurately modeling the levitation force mechanism of a spherical magnetic levitation flywheel motor. Combining the Maxwell tensor method and the equivalent magnetic circuit method, an accurate levitation force mechanism model for the spherical magnetic levitation flywheel motor is derived. Comparative analysis using finite element simulations shows that the error of the model constructed using the proposed method is reduced from 20% to 5%. This method provides an effective method, model, and approach for accurately characterizing the motor's dynamic characteristics under conditions of rotor eccentricity, magnetic circuit saturation, and diffuse magnetic leakage, laying a theoretical foundation for the analysis, design, and operational control of spherical magnetic levitation flywheel motors.

[0137] The foregoing is merely an embodiment of the present invention and is not intended to limit the present invention. It will be apparent to those skilled in the art that various modifications and variations of the present invention are possible. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention are intended to be included within the scope of the claims of the present invention.

Claims

1. A method for modeling the precise levitation force mechanism of a spherical magnetic levitation flywheel motor, characterized by: The modeling method comprises the following steps: Step 1: Obtain the magnetic field distribution of the spherical magnetic levitation flywheel motor by finite element analysis, and obtain the equivalent magnetic circuits of the radial and axial directions of the spherical magnetic levitation flywheel motor's suspension poles based on the equivalent magnetic circuit method; Step 2: Use finite element analysis to analyze the diffusion flux of the suspended pole of the spherical magnetic levitation flywheel motor. During one rotor motion cycle, based on the spherical magnetic field distribution law and the magnetic flux path, the rotor position angle range is divided into two parts. The magnetic field segmentation method is used to divide the entire air gap magnetic field into several flux conduits. The magnetic permeance of each flux conduit is calculated separately, and then the magnetic flux conduits are combined to calculate the entire air gap magnetic permeance. Step 3: Use step 2 to establish an equivalent magnetic circuit of a spherical magnetic levitation flywheel motor taking into account the diffuse magnetic flux. By analyzing and calculating the change in the air gap length between the rotor teeth and the suspended pole shoe teeth due to rotor eccentricity, the air gap magnetic permeability considering the rotor eccentricity is obtained. In addition, considering the need for magnetic circuit saturation, the air gap main magnetic flux density is corrected to obtain the corrected air gap main magnetic flux density. Step 4: Based on the equivalent magnetic circuit model of the air gap permeability obtained in step 2, combined with the air gap length and air gap main flux density obtained in step 3, analyze the different rotor position angle ranges, calculate the air gap permeability of each part, and establish the spherical magnetic levitation flywheel motor suspension force model; wherein: In step 3, the effect of rotor eccentricity on air gap permeability is considered. When the rotor is eccentric b mm in the negative direction of the y-axis, the air gap length relationship is: Further deduction: δ1≈δ0+bsinθ m When the rotor is offset a mm in the positive direction of the x-axis, the air gap length relationship is: Further deduction: δ1≈δ0+acosθ m When the rotor is offset a mm in the positive direction of the x-axis and b mm in the negative direction of the y-axis, the relationship between the air gap lengths is: Further deduction: δ1≈δ0+acosθ m -bsinθ m Among them, δ0 is the average air gap length when the rotor is not eccentric, δ1 is the air gap length between the stator and rotor poles after eccentricity, and θ m is the mechanical position angle, r is the radius of the sphere where the stator is located; The main magnetic flux density of the air gap in step 3 is: The air gap diffusion flux density is: Combined with the finite element analysis results, the air gap magnetic flux correction formula is designed as follows: Where N is the number of turns of the suspension winding, I0 is the current provided by the permanent magnet, I c is the control current of the suspension winding; B nl is the air gap flux density after correction; B c B is the critical value of the magnetic density of ferromagnetic materials tending to saturation; l According to formula B m 、B q The air gap magnetic flux density is obtained; k is the correction coefficient.

2. The method for accurately modeling the levitation force mechanism of a spherical magnetic levitation flywheel motor according to claim 1, characterized in that: In step 2, the period of the rotor movement is 30°.

3. The method for accurately modeling the levitation force mechanism of a spherical magnetic levitation flywheel motor according to claim 2, characterized in that: One cycle of rotor motion is divided into two parts: -7.5°<θ<7.5° and -15°<θ<-7.5° and 7.5°<θ<15°.

4. The method for accurately modeling the levitation force mechanism of a spherical magnetic levitation flywheel motor according to claim 1, characterized in that: In step 2, the shapes of the magnetic flux tubes divided by the magnetic field division method are divided into three shapes: rectangle, 1 / 4 solid cylinder and 1 / 4 hollow cylinder.

5. The method for accurately modeling the levitation force mechanism of a spherical magnetic levitation flywheel motor according to claim 4, characterized in that: The magnetic permeance calculation formulas of the three shapes of magnetic flux tubes are: The magnetic permeance of the rectangular magnetic flux tube under the magnetic pole is: Magnetic permeance of a 1 / 4 solid cylindrical flux tube: Magnetic permeance of a 1 / 4 hollow cylindrical flux tube: Where μ0 is the vacuum magnetic permeability, η is the arc length of the rotor tooth corresponding to the sphere, R is the radius of the rotor tooth corresponding to the sphere, h1 is the axial height of the spherical rotor tooth, δ av1 is the average air gap length in region 1, β is the angle between the axial rotor tooth spherical surface, A1 is the average cross-sectional area of ​​the rectangular flux tube in region 1, V is the volume of 1 / 4 solid cylindrical flux tube, A2 is the average cross-sectional area of ​​1 / 4 hollow cylindrical flux tube, δ av2 is the average magnetic flux length in region 2, δ av3 is the average magnetic flux length in region 3.

6. The method for accurately modeling the levitation force mechanism of a spherical magnetic levitation flywheel motor according to claim 1, characterized in that: In step 4, when the rotor position angle is -7.5°<θ<7.5°, the radial spherical air gap reluctance is: The magnetic flux generated by the ±x-direction branch and the ±y-direction branch is: Among them, G zax1 , G zax2 , G zay1 , G zay2 They are the total magnetic permeability of the air gap in the ±x and ±y directions of phase A, N is the number of turns of the A-phase suspension winding, I0 is the current provided by the permanent magnet, I c is the control current of the A-phase suspension winding, C r is the magnetic flux leakage coefficient of the radial part of the motor, Φ x1 , Φ x2 They are the magnetic flux generated by the ±x direction branches, Φ y1 , Φ y2 They are the magnetic flux generated by the ±x direction branches, in Weber (Wb), R xaq1 、R xaq2 is the air gap reluctance of phase A in the ±x direction, R yaq1 、R yaq2 is the air gap reluctance of phase A in the ±y direction; When the rotor position angle is -7.5°<θ<7.5°, the electromagnetic force F x for: Where α is half of the spherical angle of the suspended pole, and A represents the area of ​​the pole end surface including the diffusion flux range; When the rotor position angle is 15°<θ<-7.5° and 7.5°<θ<15°, the semicircular integral path dL is taken for the irregular diffusion flux at both ends and the electromagnetic force dF is decomposed to obtain the radial force dF. r and tangential force dF t for: The semicircular integral path is divided into 2n parts, n parts on each side of the symmetry axis. The radial and tangential components of the electromagnetic force generated on this path are decomposed into: When n approaches infinity, the radial and tangential components are: in, h is the thickness of the suspension of phase A or phase B, i is the number of parts, l 78 、l 56 is the integral path length; μ0 is the vacuum permeability; Therefore, the radial suspension force generated by the irregular magnetic flux at both ends of the positive x-axis of phase A is: The total radial force generated by this part in the x-axis direction of phase A is: Among them, F r2 The radial suspension force generated by the irregular magnetic flux at both ends of the negative x-axis of phase A is obtained similarly to that in the positive direction; B f , B f1 , B x1f1 , B x2f1 , B x1f , B x2f According to The formula obtains the corrected edge air gap flux density in the ±x direction; Therefore, when the rotor position angle is 15°<θ<-7.5° and 7.5°<θ<15°, the suspension force expression is: F=F x +F r 。

Images