Computer implementation method, system, and computer program for eddy current analysis

Abstract

This technology provides a preprocessing method that differs from the conventional approach of employing a single preprocessing technique. [Solution] The disclosed method is a computer implementation method performed by a computer. The computer implementation method includes numerically analyzing eddy currents by solving the discretized governing equations using a preprocessing algorithm with the finite element method. The preprocessing in the preprocessing algorithm is a composite preprocessing that performs different types of preprocessing on multiple parts of the governing equations.

Description

[Technical Field]

[0001] This disclosure relates to a computer implementation method, system, and computer program for eddy current analysis. [Background technology]

[0002] Patent Document 1 discloses a method for numerically analyzing eddy currents using the finite element method.

[0003] When numerically analyzing eddy currents using the finite element method, the discretized governing equations are solved using a preprocessing algorithm. Preprocessing is an operation performed to improve convergence or reduce the number of iterations, and is repeated during the iterations of the algorithm that solves the governing equations. Preprocessing is effective in improving the performance of the algorithm. Conventional preprocessing methods include, for example, incomplete Cholesky decomposition preprocessing. Conventional algorithms with preprocessing include, for example, the conjugate gradient method with incomplete Cholesky decomposition preprocessing. [Prior art documents] [Patent Documents]

[0004] [Patent Document 1] Japanese Patent Publication No. 2004-326630 [Non-patent literature]

[0005] [Non-Patent Document 1] Hiroo Katagiri, Kazuki Senba, and Takashi Yamada, "A Study on the Degradation of Convergence of Linear Iterative Solution Methods in Finite Element Analysis of Electromagnetic Fields Coupled with Electrical Circuits," Proceedings of the Joint Workshop on Static and Rotating Machinery, SA-21-011, RM-21-011, March 4, 2021. [Non-Patent Document 2] H. Igarashi and N. Yamamoto, “Effect of Preconditioning in Edge-Based Finite-Element Method,” IEEE Trans. Magn., vol.44, no.6, pp.942-945, 2008. [Non-Patent Document 3] Shingo Hiruma, Ken Mifune, and Tetsuji Matsuo, "Mathematical Study on the Evaluation of Condition Numbers in Darwinian Approximation of Maxwell's Equations," Proceedings of the Joint Workshop of the Institute of Electrical Engineers of Japan on Static and Rotating Machinery, SA-24-003, RM-24-003, March 2024. [Non-Patent Document 4] R. Hiptmair, “Multigrid method for Maxwell's equations,” SIAM J. Numer. Anal., vol.36, pp. 204-225, Dec. 1998. [Non-Patent Document 5] Takeshi Iwashita, Ken Mifune, and Masaaki Shimazaki, "Basic Concepts of Explicit and Implicit Error Correction Methods," Proceedings of the Research Institute for Mathematical Sciences, Vol. 1638, 18-27, 2009. [Overview of the Initiative]

[0006] In addition to the incomplete Cholesky decomposition mentioned above, several preprocessing methods have been conventionally known, including algebraic multigrid methods and domain decomposition methods. However, each of these methods has its own advantages and disadvantages, and in some cases, adopting only one of these methods as a preprocessing step may not sufficiently improve the algorithm's performance.

[0007] Therefore, there is a need for a technology to perform preprocessing using an approach different from the conventional approach of employing a single method as the preprocessing step. Accordingly, this disclosure aims to provide a technology for performing preprocessing using an approach different from the conventional approach.

[0008] The technology of this disclosure includes numerically analyzing eddy currents by solving the discretized governing equations using a preprocessing algorithm based on the finite element method. The preprocessing in the preprocessing algorithm may be a composite preprocessing that performs different types of preprocessing on multiple parts of the governing equations. In this case, different types of preprocessing can be performed depending on the part of the governing equations by the composite preprocessing that performs different types of preprocessing.

[0009] Further details will be described in the embodiments below. [Brief explanation of the drawing]

[0010] [Figure 1] Figure 1 is a diagram showing the configuration of the system according to the embodiment. [Figure 2] Figure 2 shows the governing equations. [Figure 3] Figure 3 illustrates the finite element formulation. [Figure 4] Figure 4 is a flowchart showing numerical analysis using a conventional preprocessing algorithm. [Figure 5] Figure 5(A) shows the coil wire model, Figure 5(B) shows the governing equations, and Figure 5(C) shows the estimation formula for the non-zero minimum singularity. [Figure 6] Figure 6(A) shows the basic model, and Figure 6(B) shows a model in which the length of the basic model is an integer multiple. [Figure 7] Figure 7 is a graph showing the length dependence of the non-zero minimum singularity. [Figure 8] Figure 8 explains why the Laplace equation part is a factor in the deterioration of convergence. [Figure 9] Figure 9(A) shows the composite preprocessing matrix, and Figure 9(B) is a flowchart of the numerical analysis using the algorithm with composite preprocessing. [Figure 10] Figure 10(A) shows the unknowns for the 5.5-turn model, Figure 10(B) shows the experimental results for the 5.5-turn model, and Figure 10(C) shows the experimental results for the 21.5-turn model. [Figure 11]Figure 11(A) shows the convergence characteristics for the 21.5-turn model, and Figure 11(B) is a magnified view of a portion of Figure 11(A). [Figure 12] Figure 12 illustrates the explicit error correction method. [Figure 13] Figure 13 illustrates multiple independent conductors. [Figure 14] Figure 14 illustrates parallelization. [Figure 15] Figure 15 is a conceptual diagram of hierarchical parallelization. [Figure 16] Figure 16 shows a quarter-model with 18 turns of winding. [Figure 17] Figure 17 shows the verification results of parallelization in the 18-turn winding 1 / 4 model. [Modes for carrying out the invention]

[0011] <1. Overview of computer implementation methods, systems, and computer programs for eddy current analysis>

[0012] (1) The method according to the embodiment is a computer implementation method executed by a computer, which may include numerical analysis of eddy currents by solving the discretized governing equations with a preprocessing algorithm using the finite element method. The preprocessing in the preprocessing algorithm may be a composite preprocessing that performs different types of preprocessing on multiple parts of the governing equations. In this case, different types of preprocessing can be performed depending on the part of the governing equations by a composite preprocessing that performs different types of preprocessing. Therefore, each preprocessing can function complementaryly to one another, which is preferable to conventional approaches that employ a single method as preprocessing.

[0013] (2) The plurality of parts may include a first part containing a discretized form of Ampere's law and a second part containing a discretized form of the current continuity equation. The combined preprocessing includes a first preprocessing and a second preprocessing of a different type from the first preprocessing, wherein the first preprocessing may be a preprocessing for the first part which is a discretized form of Ampere's law. The second preprocessing may be a preprocessing for the second part which is a discretized form of the current continuity equation. In this case, convergence can be improved.

[0014] (3) Preferably, the second pretreatment is a pretreatment that can make the lower bound of the minimum singularity larger than when the first pretreatment is used.

[0015] (4) The first pretreatment may be an incomplete Cholesky decomposition pretreatment.

[0016] (5) The second preprocessing step may be based on an algebraic multigrid method or a domain decomposition method.

[0017] (6) The preprocessing algorithm may be a preprocessing conjugate gradient method.

[0018] (7) The eddy current may be an eddy current inside the coil winding.

[0019] (8) The computer implementation method may further include providing a user interface for setting parameters for each preprocessing that are set to execute each preprocessing included in the composite preprocessing, and accepting the setting of the parameters for each preprocessing through the user interface. Each preprocessing included in the composite preprocessing may be executed based on the parameters set for each preprocessing.

[0020] (9) The eddy current analysis may include a plurality of independent conductor regions. The plurality of parts may include a part that includes a discretized equation of the current continuity equation. The combined preprocessing may include preprocessing for the part that includes the discretized equation of the current continuity equation. The preprocessing for the part that includes the discretized equation of the current continuity equation may have a plurality of preprocessing steps corresponding to the plurality of independent conductor regions. The plurality of preprocessing steps corresponding to the plurality of independent conductor regions may be processed in parallel. Processing can be sped up by parallel processing.

[0021] (10) The combined preprocessing may include a first preprocessing and a second preprocessing of a different type from the first preprocessing. The first preprocessing and the second preprocessing may be processed in parallel. Parallel processing can speed up the processing.

[0022] (11) The eddy current analysis may include a plurality of independent conductor regions. The plurality of parts may include a first part including a discretized equation of Ampere's law, and a second part including a discretized equation of the current continuity equation. The combined preprocessing may include a first preprocessing for the first part and a second preprocessing for the second part. The first and second preprocessing may be processed in parallel. The second preprocessing may have a plurality of preprocessing corresponding to the plurality of independent conductor regions. The plurality of preprocessing corresponding to the plurality of independent conductor regions may be processed in parallel. Processing can be sped up by parallel processing.

[0023] (12) The system according to the embodiment may be a system for numerically analyzing eddy currents by solving the discretized governing equations using the finite element method with a preprocessing algorithm. The preprocessing in the preprocessing algorithm may be a composite preprocessing that performs different types of preprocessing on multiple parts of the governing equations.

[0024] (13) The computer program according to the embodiment may be a computer program that causes a computer to perform processing for numerical analysis of eddy currents by solving the discretized governing equations using the finite element method with a preprocessing algorithm. The preprocessing in the preprocessing algorithm may be a composite preprocessing that performs different types of preprocessing on multiple parts of the governing equations.

[0025] The computer program may be stored in a computer-readable, non-temporary storage medium.

[0026] <2. Examples of computer implementation methods, systems, and computer programs for eddy current analysis>

[0027] The embodiments will be described in more detail below with reference to the drawings.

[0028] Figure 1 shows a system 100 according to an embodiment. System 100 performs numerical analysis using the finite element method. The processing performed by system 100 may include numerical analysis of eddy currents occurring in the object under analysis. Eddy current analysis is performed, for example, for the product design of drive motors used in hybrid vehicles and electric vehicles. Eddy current analysis is also performed in the product design of stationary equipment such as transformers.

[0029] System 100 may be configured with one or more computers. As an example, System 100 shown in Figure 1 is configured with a computer comprising a processor 10 and memory 20 connected to the processor 10. The processor 10 is, for example, a CPU, GPU, or other type of processor. The memory 20 may comprise, for example, primary and secondary storage devices. The primary storage device is, for example, RAM. The secondary storage device is, for example, a hard disk drive (HDD) or a solid-state drive (SSD).

[0030] Memory 20 may contain a computer program 21 executed by the processor 10. The processor 10 reads and executes the computer program 21 stored in memory 20. The computer program 21 has program code that indicates instructions to the computer to realize some or all of the functions of the system 100 according to this embodiment. Memory 20 can also store parameters 22, which will be described later. Memory 20 can also store data to be analyzed 23. The system 100 performs numerical analysis of the data to be analyzed 23 using the finite element method.

[0031] The computer comprising System 100 may be equipped with an input / output interface 30. An input device 50 and an output device 60 may be connected to the input / output interface 30. The input device 50 is a device for input, such as a keyboard or mouse. System 100 can acquire data from an external source via the input device 50. The output device 60 is a device for output, such as a display. System 100 can provide a user interface screen to the user of System 100 via the output device 60. The input device 50 and the output device 60 may be integrated into a single unit, such as a touch panel display.

[0032] The computer comprising System 100 may be equipped with a communication interface 40. The communication interface 40 is an interface for communicating with external devices via a network. System 100 can acquire data from external sources via the communication interface 40. Furthermore, System 100 can provide a user interface screen to the display (output device 60) of a terminal used by a user of System 100 via the communication interface 40.

[0033] <2.1 Governing equation>

[0034] Figure 2 illustrates the formulation of the governing equations for the eddy current problem analyzed by System 100. Figure 2 shows both the AV method (or A-φ method) formulation, which utilizes the magnetic vector potential and electric scalar potential, and the A method formulation, which generally involves fewer unknowns than the AV method. Note that in Figure 1 and other figures, as well as in the following explanation, symbols not otherwise specified are used according to the usual notation in the field of electromagnetism (for the usual notation usage, see, for example, Non-Patent Documents 1-5).

[0035] Discretization is necessary to numerically solve the governing equations. There are various formulations for discretization, but for electromagnetic devices driven at low frequencies, the formulation using the finite element method is the most widely used. The finite element method is a technique that discretizes space by dividing it into a computational grid called a finite element and assigning unknowns to the edges or nodes. In the case of electromagnetic fields, the edge element finite element method is used to represent the continuity of physical quantities at the material boundaries (continuity of the tangential component of the magnetic field, continuity of the normal component of the magnetic flux density).

[0036] <2.2 Finite Element Formulation and Algorithms with Preprocessing>

[0037] Figure 3 illustrates the finite element formulation. Figure 3 shows the governing equations discretized for the finite element method. These governing equations are a system of equations obtained by discretizing Ampere's law and the equation for current continuity. As explained in Figure 3, the governing equations discretized for the finite element method are solved using an iterative algorithm. A pre-processed algorithm is used as the iterative method. An example of a pre-processed algorithm is the pre-processed conjugate gradient method.

[0038] As mentioned above, method A generally has fewer unknowns than method AV. However, method A is disadvantageous from the perspective of computation time because the number of iterations required for convergence increases when using an iterative solution method. Therefore, in industry, it is common to use method AV, which can accelerate convergence. By applying preprocessing to method AV, the lower bound of the non-zero minimum singular value included in the coefficient matrix becomes a finite value that does not depend on frequency, thereby improving convergence.

[0039] When using the conjugate gradient method as an iterative solution algorithm, the choice of preprocessing for the matrix is ​​crucial for stable convergence and reduction of iterations. Preprocessing refers to operations that improve the condition number of the equation and are performed repeatedly during the iterations of the iterative solution algorithm, such as the conjugate gradient method. The condition number can be expressed, for example, as the ratio of the largest singular value to the smallest singular value of the matrix corresponding to the equation.

[0040] Figure 4 shows an iterative solution algorithm with preprocessing for numerically solving the governing equations. As shown in Figure 4, in the iterative solution algorithm with preprocessing, preprocessing using the preprocessing matrix M is repeatedly performed during the iterations in step S20, which is iterated through to solve the governing equations.

[0041] The preprocessing matrix M used in preprocessing is obtained by approximating the coefficient matrix A (A ≈ M) and taking the inverse operation (M -1 It is desirable to use a matrix with a low computational cost. For matrix A, choosing the preprocessing matrix M=A is the best preprocessing method, as it results in only one iteration until convergence.

[0042] However, since it is practically impossible to calculate the inverse matrix, we usually consider an approximate inverse. The simplest preprocessing is diagonal preprocessing. This involves using the following preprocessing matrix: M=diag(A) We choose this. Calculating the inverse of a diagonal matrix is ​​easy, so it is often used in educational contexts. However, diagonal preprocessing has a weak preprocessing effect, and convergence improvement cannot often be expected for difficult problems. As a computationally inexpensive and versatile preprocessing method, there is the incomplete Cholesky (IC) preprocessing. In IC decomposition, the coefficient matrix A is A=LDL T +R This is how it is decomposed. Here, L is a lower triangular matrix. In the usual Cholesky decomposition, L is a dense matrix and R=0, but in the case of a large sparse matrix, the dense matrix cannot be preserved. Therefore, in order to reduce the computational cost of L, we reduce the computational cost while ensuring approximation accuracy by mimicking the non-zero patterns of the sparse matrix A. In this case, M=LDL T This is the pretreatment for incomplete Cholesky decomposition.

[0043] More computationally intensive and robust preprocessing methods include the Domain decomposition method (DDM) and the Algebraic multi-grid method (AMG).

[0044] <2.3 Composite Preprocessing>

[0045] The system 100 according to this embodiment performs a composite preprocessing as a preprocessing step. Composite preprocessing is a preprocessing step that combines different types of preprocessing. In composite preprocessing, the type of preprocessing applied to one part of the governing equation differs from the type of preprocessing applied to other parts. Therefore, in composite preprocessing, instead of applying a single preprocessing step to the entire governing equation, different types of preprocessing can be performed depending on the part of the governing equation. As a result, it becomes possible to apply appropriate preprocessing steps according to the part of the governing equation.

[0046] While several preprocessing methods have been known for some time, adopting only one of these methods may not adequately improve algorithmic performance in some cases. However, in combined preprocessing, different types of preprocessing can complement each other, making it preferable to the conventional approach of adopting only one preprocessing method.

[0047] The types and combinations of preprocessing steps that constitute the composite preprocessing are not particularly limited. The types and combinations of preprocessing steps that constitute the composite preprocessing can be appropriately selected depending on the individual problem in which the preprocessing may be applied. For example, the preprocessing steps that constitute the composite preprocessing can be selected from the group consisting of incomplete Cholesky decomposition, preprocessing based on domain decomposition, and algebraic multigrid methods (AMG).

[0048] The system 100 can provide the user with a user interface for setting parameters for each preprocessing step included in the composite preprocessing. For example, the system 100 can display a user interface screen for setting parameters 22 on the output device 60, which is a display. The user interface screen is a screen on which the user performs operations to set parameters 22. In the system 100 according to this embodiment, each of the multiple preprocessing steps included in the composite preprocessing is executed, so it is preferable to be able to set parameters 22 for each preprocessing step.

[0049] System 100 receives user input from the user via the input device 50 to operate the user interface screen, accepts the setting of parameters 22 for each preprocessing step, and stores the parameters 22 in memory 20. System 100 executes each preprocessing step included in the composite preprocessing step based on the parameters 22 set for each preprocessing step.

[0050] If the preprocessing is an incomplete Cholesky decomposition, system 100 can accept the following parameter 22 settings, for example. • Acceleration factor: 1.05

[0051] If the preprocessing is an algebraic multigrid method (AMG), system 100 can accept the following parameter 22 settings, for example. • Cycle type: V or W • Number of pre-relaxation sessions (npre): 3 • Number of post-easing easing sessions (npost): 3 • Maximum number of levels (max_level): 5 • Number of cycles (pre_cycles): 2 ·Smoother: Gauss-Siedel, Jacobi, IC, etc.

[0052] If the preprocessing method is domain decomposition, the system can accept the following parameter 22 settings, for example. ·Number of divisions: N • Overlap area setting (how much overlap there should be): M • Pretreatment: Schwarz pretreatment or 2-level Schwarz method • Local solver: Direct method (LU decomposition) or iterative method (CG method or GMRES) • Coarse solver: Direct method (LU decomposition) or iterative method (CG method or GMRES) • Boundary condition handling: Dirichlet boundary, Neumann boundary, or impedance boundary

[0053] The system 100 may also be able to accept the user's specification of the type of preprocessing that constitutes the composite preprocessing. In this case, the user can freely specify the combination of preprocessing that constitutes the composite preprocessing.

[0054] In the following, as an example, we will describe a case in which a complex preprocessing is used to address the problem reported in Non-Patent Document 1 (which requires a very long time to obtain a solution (resulting in poor convergence)) in the numerical analysis of an eddy current problem using the finite element method.

[0055] System 100, which performs complex preprocessing to address the above problems, can accelerate eddy current analysis. Therefore, for example, it can quickly calculate eddy current losses occurring in the analysis target, such as a design model. Furthermore, it enables many trial-and-error iterations within a short design period. Moreover, this technology can be used to perform calculations even for problems of a size that were previously difficult to calculate, thereby realizing higher-performance electromagnetic equipment designs.

[0056] <2.4 Problems in Non-Patent Document 1>

[0057] This document outlines the problems described in Non-Patent Document 1. Non-Patent Document 1 reports that the convergence of the linear iterative solution method deteriorates drastically in the analysis of eddy current losses occurring in coil wires. The analysis of eddy current losses occurring in coil wires is important for estimating AC copper losses (AC losses) caused by high-frequency currents in PWM control used in rotating machinery and the like. Non-Patent Document 1 reports that convergence deterioration occurs in various models, and the composite preprocessing according to the embodiment can, in principle, improve the convergence characteristics for all of them.

[0058] For the coil wire model shown in Figure 5(A), we assume a voltage source at both ends of the coil wire and apply a forced current. The governing equation is expressed as equation (1) in Figure 5(B).

[0059] For the model in Figure 5(A), the convergence of the linear iterative solution method deteriorates significantly, and a solution cannot be obtained when the number of turns is large. This is a problem that occurs in a great many cases in the analysis of actual machines, and a solution is desired. Numerical experiments have revealed the following experimental facts (1) to (4).

[0060] (1) Increasing the number of turns increases the number of iterations in roughly proportion to the number of turns. (2) Increasing the mesh size also increases the number of iterations. (3) Convergence deteriorates in both the time domain and the frequency domain. (4) In the static magnetic field limit ω->0, there is a gap because it does not match the number of iterations of the static magnetic field analysis.

[0061] Comparing the above experimental facts (1) to (4) with the estimation formula for the non-zero minimum singular value in equation (2) of Figure 5(C), the following can be inferred. Note that equation (2) is the estimation formula for the non-zero minimum singular value for the equation obtained by applying block diagonal preprocessing to the expanded equation of the pseudo-singular equation newly obtained in Non-Patent Document 3.

[0062] (Speculation) When the eddy current problem is formulated using the AV method and the discrete equations are considered, the matrix G corresponds to the coefficient matrix of Laplace's equation. T The smallest singular value of [σ]G corresponds to matrix C, which is the coefficient matrix of the equation for the static magnetic field. T It is smaller than the minimum singularity of [ν]C, and the former depends on the mesh size and the number of turns.

[0063] <2.5 Preliminary Experiment>

[0064] Based on the above inferences, preliminary experiments were conducted to calculate the minimum singular values ​​of the coefficient matrices for the Laplace equation, the equation for static magnetic fields (MS), and the eddy current problem (MQS) for the model shown in Figure 6, and to calculate their dependence on the length of the conductor (coil). Figure 6(A) shows the basic model of the conductor, and Figure 6(B) shows a model in which the length of the basic model is an integer multiple.

[0065] Figure 7 shows the results of a preliminary experiment, illustrating the length dependence of the non-zero minimum singularity. From these results, the Laplace equation is expressed as follows with respect to length: L -2 It can be seen that there is a dependence on length. Where the singular values ​​of Laplace's equation are large, the singular values ​​of the eddy current problem (MQS) are flat, which is the same characteristic as the singular values ​​of the equation for a static magnetic field (MS). Furthermore, where the singular values ​​of Laplace's equation are small, it can be seen that the singular values ​​of the eddy current problem (MQS) also have a dependence on length.

[0066] Thus, it can be seen that the dependence of the singular values of the eddy current problem (MQS) is determined according to the magnitudes of the singular values of the magnetostatic field equation (MS) and the Laplace equation. This result does not conflict with the prediction of the singular values expected from the estimation formula of Equation (2) in Fig. 5(C).

[0067] The experimental facts shown by the preliminary experiment mean the following.

[0068] (1) The dependence on the size of the region of the magnetostatic field equation (MS) is not strong, and it does not become a bad condition even in a relatively wide region. (2) Since the Laplace equation has a strong dependence on the size of the region and the mesh size, when making a long conductor or a fine mesh, the condition number extremely easily deteriorates, which brings about the ill-conditioning of the equation of the eddy current problem (MQS).

[0069] As also explained in Fig. 8, the minimum singular value of the Laplace equation depends on the mesh and also on the length L (the size of the analysis region). Also, the minimum eigenvalue of the Laplace equation is proportional to L -2 . In addition, due to the combined mesh dependence of the Laplace equation, the deterioration of convergence in eddy current analysis is presumed to be caused by the small singular values existing in the matrix G T [σ]G. Therefore, for the improvement of convergence, a prescription for the matrix G T [σ]G is effective.

[0070] Furthermore, the problem of degraded convergence described above can be attributed to the equation for current continuity. Therefore, by applying different pretreatments to Ampere's law and the equation for current continuity in the governing equations, the condition number of the coefficient matrix can be improved, thereby improving the convergence characteristics. Conventionally, a single pretreatment (such as incomplete Cholesky decomposition) was used for the entire matrix, which failed to effectively utilize the inherent nature of the problem in the equations, leading to degraded convergence. While AMG can also be applied to the entire coefficient matrix to improve convergence, it degrades parallel performance, making it difficult to reduce computation time. According to the technology of this embodiment, it is possible to improve the efficiency of solving while maintaining parallel performance, and as a result, computation time can be reduced.

[0071] <2.6 Use of combined preprocessing (block diagonal preprocessing) to improve convergence>

[0072] Preliminary experiments show that modifying the part of the current continuity equation in the governing equation is effective in improving convergence. Note that, below, "part" may be referred to as "block." "Block" may be used to refer to a portion of a matrix.

[0073] Matrix G included in the coefficient matrix of the eddy current problem T To improve the condition number of [σ]G, the algebraic multigrid method (AMG) is effective. On the other hand, matrix C T Algebraic multigrid methods for [ν]C are highly dependent on the convergence characteristics of the problem, and in some cases, convergence may be slow.

[0074] Therefore, as a preprocessing step, matrix C T Using incomplete Cholesky decomposition pretreatment on the [ν]C block, the matrix G T The convergence can be dramatically improved by applying a composite preconditioner, such as the algebraic multigrid method (AMG), to the [σ]G block. This is also shown in equation (2) of Figure 5(C), and by performing different preconditioners on each part to which it is applied, the lower bound for the non-zero minimum singularity can be increased.

[0075] The equation for current continuity is dependent on the equation for Ampere's law, and these simultaneous equations are redundant. However, by explicitly solving the equation for current continuity, different pretreatments can be applied. Furthermore, in composite pretreatment, convergence can be improved with a simple block-diagonal pretreatment. Here, the composite pretreatment matrix for composite pretreatment is defined as a block-diagonal matrix as shown in equation (3) of Figure 9(A).

[0076] As shown in equation (3), the composite preprocessing matrix is ​​M as a block of matrices. A and M φ It has. Therefore, the composite pretreatment is M A First pretreatment by and M φ This includes a second pretreatment by the following: In formula (3), the first pretreatment M A C T This is a pretreatment matrix for [ν]C+jω[σ], and is the second pretreatment M. φ is 1 / jωG T This is the pretreatment matrix for [σ]G.

[0077] First pretreatment M A For example, an incomplete Cholesky decomposition pretreatment can be employed. Second pretreatment M φ For example, preprocessing methods based on algebraic multigrid (AMG) or domain decomposition can be used. First preprocessing M A and second pretreatment M φ Each pretreatment is free to choose from, but each pretreatment M A M φ It is preferable to select the one that is effective depending on the part to which it is applied. First pretreatment M A and second pretreatment M φVarious combinations are possible, such as "incomplete Cholesky decomposition + AMG" and "incomplete Cholesky decomposition + domain decomposition method." A theoretically interesting point is that the convergence of the current continuity equation, which is introduced to improve the convergence of the formulation by method A, can be further improved by applying strong pretreatment such as algebraic multigrid methods.

[0078] By increasing the lower bound of the minimum singularity, the condition number can be improved and convergence can be improved, thus the second preprocessing M φ Therefore, a method that can significantly increase the lower bound of the minimum singularity is preferable. A When the type is appropriately determined, the second pretreatment M φ The first pretreatment M was determined. A It is preferable to use a pretreatment that can increase the lower bound of the minimum singularity compared to when using the first pretreatment M. A However, in the case of incomplete Cholesky decomposition, the second pretreatment M φ By choosing either the algebraic multigrid method or the domain decomposition method, the condition number can be improved by increasing the lower bound of the minimum singularity. As a result, convergence can be dramatically improved.

[0079] Figure 9(B) shows an iterative solution algorithm with compound preprocessing. As shown in Figure 9(B), in the iterative solution algorithm with compound preprocessing, compound preprocessing using the compound preprocessing matrix M (equation (3)) is repeatedly performed during the iterations in step S90, which is iterated to solve the governing equations.

[0080] As mentioned above, the governing equations to which the combined pretreatment is applied are a system of equations obtained by discretizing Ampere's law and the equation for current continuity. In step S90, of the governing equations, the first pretreatment M A The part to which it applies is the first part, which includes the discretized form of Ampere's Law, and the second pretreatment M φThe part to which this applies is the second part, which includes the discretized equation for current continuity. By using different pretreatments on the system of equations obtained by discretizing Ampere's law and the equation for current continuity, the efficiency of linear iterative solving can be improved.

[0081] <2.7 Numerical experiments on combined preprocessing>

[0082] A composite preprocessing step defined by a block diagonal matrix (equation (3)) was applied to the 5.5-turn model in Figure 5(A) (Example). The unknowns of the 5.5-turn model are shown in Figure 10(A). First preprocessing step M A For this, incomplete Cholesky decomposition (IC) is used, and the second pretreatment M φ The algebraic multigrid method (AMG) was used. Numerical calculations were performed for both the V-cycle (Example 1A) and the W-cycle (Example 1B) as parameters for the algebraic multigrid method (AMG). In both Examples 1A and 1B, the Gauss-Seidel smoother and aggregation courseing were used as other parameters for the AMG, and other parameters were also used as deemed appropriate.

[0083] Figure 10(B) shows the calculation results when compound preprocessing (CP) is performed on the 5.5-turn model (Examples 1A and 1B). As a comparative example, Figure 10(B) also shows the calculation results when single preprocessing is performed using incomplete Cholesky decomposition (IC). As is clear from Figure 10(B), in Examples 1A and 1B, both the number of iterations and the calculation time were reduced to about 1 / 10 compared to the comparative example. Matrix G by AMG T The improvement in the number of conditions in [σ]G, the reduction in the number of iterations, and the significantly smaller number of nodes compared to the number of edges, which reduces the computational load of AMG, are thought to have contributed to the reduction in computation time.

[0084] Figure 10(C) shows the calculation results for Examples 2A, 2B, and 2C and the comparative example for the more complex 21.5-turn model. In Examples 2A and 2B, the first pretreatment M AAs a second pretreatment, incomplete Cholesky decomposition (IC) is used. φ The Algebraic Multigrid Method (AMG) was used as the parameter. Example 2A, like Example 1A, adopted the V cycle as the parameter for AMG, and Example 2B, like Example 1B, adopted the W cycle as the parameter for AMG. In Example 2C, the first pretreatment M A As a second pretreatment, incomplete Cholesky decomposition (IC) is used. φ Domain Decomposition Method (DDM) was used as the method. In the 21.5-turn model, in Examples 2A, 2B, and 2C, a significant reduction in both the number of iterations and computation time was achieved compared to the comparative example. Thus, by performing a combined preprocessing using AMG or DDM, the matrix G T The number of conditions in [σ]G improves, and the number of iterations decreases.

[0085] Figure 11(A) shows the convergence characteristics of Examples 2A, 2B, and 2C and the comparative example. Figure 11(B) is a magnified view of the convergence characteristics of Examples 2A, 2B, and 2C shown in Figure 11(A). It can be seen that Examples 2A, 2B, and 2C all converge much faster than the comparative example.

[0086] <2.8 Supplement>

[0087] It should be noted that the expression as a composite preprocessing method is not limited to equation (3) shown in Figure 9(A). For example, it is possible to achieve the same convergence characteristics as when using composite preprocessing by equation (3) within the framework of method A. This can be done by utilizing a framework called explicit error correction (see Non-Patent Literature 5).

[0088] Explicit error correction methods describe the algebraic relationship between linear iterative solution methods derived from the equivalence of Hybrid Smoother and AV methods; AV methods can be considered implicit error correction methods, while Hybrid Smoother is an explicit error correction method. Using this framework, a formulation of explicit error correction methods can be created to achieve convergence equivalent to that of composite preprocessing.

[0089] As explained in detail in Figure 12, even in the explicit error correction method, equivalent convergence characteristics can be obtained as when using the block diagonal matrix shown in equation (3) as the composite preprocessing matrix.

[0090] In any case, the essential point from the perspective of improving convergence is that the efficiency of linear iterative solving can be improved for the eddy current problem in windings by using different pretreatments on the system of equations obtained by discretizing Ampere's law and the equation for current continuity. In particular, when the length of the winding is long, the convergence characteristics can be improved by using a strong pretreatment on the equation for current continuity.

[0091] <2.9 Parallelization of preprocessing>

[0092] The following describes methods for parallelizing preprocessing. Here, as an example, for an analysis model (the object of analysis) that includes multiple independent conductor regions, the preprocessing is parallelized based on the independence of the physical regions (parallelization considering the independence of the conductor regions). Furthermore, the parallelization may be the parallelization of the first preprocessing and the second preprocessing included in the composite preprocessing (parallelization of each preprocessing that constitutes the composite preprocessing). In addition, the parallelization may be parallelization in the computational domain of the application of the preprocessing. The parallelization may be a combination of two or more of these parallelizations. Below, as an example, it will be shown that a hierarchical parallelization method can be constructed by combining parallelizations.

[0093] Generally, electromagnetic devices form circuits containing numerous series and parallel connections, requiring the modeling of many winding regions. In such cases, the analytical model (computational model) often includes multiple independent conductor regions. In this case, parallelization can be efficiently performed using the method described below. Here, multiple independent conductor regions may be multiple conductor regions insulated so that no current flows between them, but it is sufficient if the different conductor regions, like the multiple windings in the computational model, form independent circuits in terms of the circuit equations.

[0094] Figure 13 shows multiple independent conductor regions Ω. ciThis shows Ampere's law and the equation for current continuity in the analytical region Ω having multiple independent conductor regions Ω. ci These are, for example, multiple winding regions in electromagnetic equipment. Equation (4) in Figure 13 shows Ampere's law, and equations (5) to (7) in Figure 13 show Ω c1 From Omega cm The equations for current continuity for each conductor region up to Ω are shown. As shown in equations (5) to (7), the equations for current continuity are given for multiple independent conductor regions Ω ci Each is given. Equations (8) to (11) in Figure 13 show examples of boundary conditions given to the ends of the boundary. Note that Γ in equation (10) in,i This is the conductor region Ω ci This is the current inflow boundary in equation (11), and Γ out,i This is the conductor region Ω ci This is the current outflow boundary in that region.

[0095] Furthermore, Figure 13 also illustrates the aforementioned composite pretreatment, and equation (12) in Figure 13 is the same as equation (1) in Figure 5. Also, equation (13) in Figure 13 is the same as equation (3) in Figure 9, and represents the composite pretreatment M. The composite pretreatment M is the first pretreatment M A and second pretreatment M φ It has the following characteristics. As shown in equation (13), the composite preprocessing matrix M is represented as a block diagonal matrix.

[0096] Figure 14 shows the second pretreatment M. φ This shows the parallelization of the process and the entire composite preprocessing M. Focusing on the matrix structure of equation (12) in Figure 13, since the conductor region to be analyzed is composed of m independent regions, G in equation (12) T [σ]G has the structure of equation (14) in Figure 14. Therefore, the second preprocessing matrix M in the composite preprocessing M φ As a result, the matrix shown in equation (15) in Figure 14 is obtained. As shown in equation (15), the second preprocessing M φ This is a multiple conductor region Ω ci Multiple preprocessing steps M corresponding to i (Multiple partial preprocessing steps M) i ) has. Thus, the second preprocessing matrix M φThe self is also represented as a block diagonal matrix.

[0097] Equation (16) in FIG. 14 shows a composite preprocessing similar to step S90 in FIG. 9(B). M in equation (16) A -1 r1 k calculation and M φ -1 r2 k calculation can be performed independently by utilizing the block diagonal structure of the composite preprocessing matrix M. That is, in equation (16), M A -1 r1 k calculation and M φ -1 r2 k calculation are processed in parallel (parallelization of each preprocessing included in the composite preprocessing). By this parallelization, the calculation efficiency of the composite preprocessing can be improved.

[0098] z2 in equation (16) k = M φ -1 r2 k is represented as in equation (17) by utilizing equation (15). M1 in equation (17) -1 r 21 k , M2 -1 r 22 k ··· M m 1 r 2m k each calculation can be performed independently by utilizing the block diagonal structure of the matrix. That is, in equation (17), M1 constituting the second preprocessing -1 r 21 k , M2 -1 r 22 k ··· M m 1 r 2m kEach calculation is processed in parallel (parallelization considering the independence of the conductor domain). This parallelization improves the computational efficiency of the second preprocessing step. Moreover, even for preprocessing methods that generally do not have high parallelization efficiency, such as AMG or direct methods, parallelization is possible for each independent physical domain without impairing the convergence improvement effect of the preprocessing.

[0099] Then, from equations (16) and (17), equation (18) is obtained. M in equation (18) A -1 r1 k ,M1 -1 r 21 k ,M2 -1 r 22 k ···M m 1 r 2m k Each of these calculations can be performed independently by utilizing the structure of the block diagonals of the matrix. That is, in equation (18), M A -1 r1 k ,M1 -1 r 21 k ,M2 -1 r 22 k ···M m 1 r 2m k Each of these calculations can be performed in parallel. According to equation (18), the first preprocessing M A and second pretreatment M φ This can be parallelized, and furthermore, the second preprocessing M φ Parallelization of each computation that constitutes the process, as well as combined parallel processing, are possible. Combined parallelization can further improve the computational efficiency of the complex preprocessing.

[0100] To compare with the case where the independence of the conducting region is not considered, let's consider the case where the direct method is used as the second preprocessing matrix. In this case, the preprocessing matrix is ​​M as the product of an upper triangular matrix and a lower triangular matrix. φ =LL TIt is written as follows. Applying this preprocessing matrix is ​​equivalent to performing forward and backward substitution on the residual vector. Forward and backward substitution must be calculated sequentially, and their parallelization efficiency is generally not high. On the other hand, when the independence of the conductor regions is considered, the preprocessing for the m conductor regions is considered independently, so the parallelization efficiency is extremely high.

[0101] Moreover, M A -1 r1 k ,M1 -1 r 21 k ,M2 -1 r 22 k ···M m 1 r 2m k Each of the calculations themselves can also be subdivided and parallelized (parallelization within the computational domain). By combining this with parallelization within the computational domain, the computation time can be further reduced.

[0102] As described above, by considering all aspects of the preprocessing of each process constituting the composite preprocessing, parallelization that takes into account the independence of the conductor domain for the second preprocessing (parallelization based on the physical domain), and parallelization in the computation domain to which the preprocessing matrix is ​​applied, it becomes possible to schedule hierarchical parallelization according to the available computing resources. Figure 15 conceptually illustrates an example of such hierarchical parallelization.

[0103] In Figure 15, M A -1 r1 k ,M1 -1 r 21 k ,M2 -1 r 22 k ···M m 1 r 2m k For each calculation, one or more computing nodes are allocated as appropriate, taking into account the computational load, etc. The multiple computing nodes to which each calculation is allocated execute the calculation in parallel. This allows M A-1 r1 k ,M1 -1 r 21 k ,M2 -1 r 22 k ···M m 1 r 2m k Each calculation is processed in parallel across multiple computing nodes. Furthermore, the calculations within each computing node are subdivided and executed in parallel by multiple cores (parallelization of the computing domain). A computing node is, for example, a processor such as a CPU or GPU. A core is, for example, an arithmetic unit (core) in a multi-core processor. Each computing node can execute multiple processes in parallel using multiple cores.

[0104] <2.9 Numerical Verification of Parallelization>

[0105] In the following, we consider a model (analytical symmetry) in which an 18-turn coil is wound around a magnetic core as the object of analysis for verifying parallelization efficiency. Here, due to the symmetry of the model, we use a 1 / 4 symmetric model as shown in Figure 16. Note that the model to be analyzed may be any symmetric model other than the 1 / 4 symmetric model. Also, the model to be analyzed does not have to be a symmetric model that is a part of the analysis model, as long as it contains multiple independent conductive regions.

[0106] The windings are assumed to be wound in series. The windings consist of a conductor surrounded by insulating material. In the 1 / 4 symmetric model shown in Figure 16, there are multiple independent conductor regions Ω. ciThe system has multiple winding regions (18 winding regions: i=1 to 18). Each winding (conductor region) is insulated from the others. Because it is a symmetric model (axisymmetric model), both longitudinal ends of each of the 18 windings are conductor cross-sections. The cross-section at one longitudinal end of each winding is the current inflow surface, and the cross-section at the other end is the current outflow surface. Ampere's law and the continuity equation are given by equations (4) to (7) in Figure 13, and the boundary conditions are given by equations (8) to (11). Here, finite element equations were created under the condition that a current of 1 [A] flows through the cross-section of each conductor (winding).

[0107] The 18 windings in Figure 16 are each physically independent, and as mentioned above, they constitute multiple independent conductor regions. In a symmetric model, the model itself does not have connection information for conductors outside the boundaries. However, each winding is in a series connection relationship when considering the areas outside the boundaries. Even though they are in a series connection relationship, they can be treated as independent depending on how the unknowns placed in the conductor cross-sections are coupled with the circuit equations. In other words, in the circuit equations representing the model in Figure 16, if each conductor region is represented as an independent circuit, and only the current value at the current inflow cross-section of each conductor is specified, then, according to the current continuity equation, the same value of current will automatically flow out at the current outflow cross-section. Therefore, if the same current value is input at each cross-section of multiple conductors, the same current will flow through all conductors as a result. For this reason, even if multiple independent conductor regions form independent circuits in the circuit equations, a series connection can be represented. Thus, when different conductor regions form independent circuits in the circuit equations, even if multiple conductor regions are connected in series, they can be treated as if multiple independent conductor regions exist.

[0108] The mesh generated for the model shown in Figure 16 had 559,976 elements, 828,300 edges, and 167,880 nodes. The discretized equations were solved using a composite preprocessing method, employing incomplete Cholesky decomposition (IC) as the first preprocessing step and the direct method (Direct) as the second preprocessing step. The specifications of the computer used for the calculations were: CPU: i9-10980XE CPU@3.00GHz, Memory: DDR4-2933, rank 1, 8GBx8=64GB.

[0109] We compared two methods: one that parallelizes 18 independent conductors and assigns one core to each conductor domain to apply the direct method (Example B: Physical Domain Parallelization), and another that applies the direct method to all conductor domains together and parallelizes forward and backward substitution (Example A: Computation Domain Parallelization). Convergence was judged using a relative residual norm of 1e-8. The results are shown in Figure 17. The calculation results show that the computation time is significantly reduced by physical domain parallelization. Furthermore, there was no change in the number of iterations depending on the parallelization method, indicating that it is possible to reduce computation time while maintaining the convergence improvement effect. In addition, the number of iterations was 630 for both Examples A and B, which is fewer than the comparative examples shown in Figures 10 and 11, indicating improved convergence compared to the comparative examples.

[0110] The present invention is not limited to the above embodiments, and various modifications are possible. [Explanation of Symbols]

[0111] 10: Processor 20: Memory 21: Computer Programs 22: Parameters 23: Data to be analyzed 30: Input / Output Interface 40: Communication Interface 50: Input Devices 60: Output device 100: System

Claims

1. A computer implementation method performed by a computer, the computer implementation method includes numerically analyzing eddy currents by solving the discretized governing equations using a preprocessing algorithm with the finite element method, The preprocessing in the aforementioned preprocessing algorithm is a composite preprocessing that performs different types of preprocessing on multiple parts of the governing equation. Computer implementation methods for eddy current analysis.

2. The aforementioned parts are, Part 1 includes the discretized form of Ampere's Law, The second part includes the discretized equation for current continuity, Includes, The aforementioned combined pretreatment includes a first pretreatment and a second pretreatment of a different type from the first pretreatment. The first pretreatment is a pretreatment of the first part, which is a discretized form of the equation of Ampere's law, The second pretreatment is a pretreatment of the second part, which is a discretized equation of the current continuity equation. The computer implementation method according to claim 1.

3. The second preprocessing step is a preprocessing step that can increase the lower bound of the minimum singularity compared to when the first preprocessing step is used. The computer implementation method according to claim 2.

4. The first pretreatment described above is an incomplete Cholesky decomposition pretreatment. The computer implementation method according to claim 2.

5. The second preprocessing step is a preprocessing step based on an algebraic multigrid method or a domain decomposition method. The computer implementation method according to claim 2.

6. The aforementioned preprocessed algorithm is a preprocessed conjugate gradient method. The computer implementation method according to claim 1.

7. The eddy currents mentioned above are eddy currents inside the coil windings. The computer implementation method according to claim 1.

8. The computer implementation method provides the user with a user interface for setting parameters for each preprocessing step included in the composite preprocessing step, The user interface accepts the setting of the parameters for each preprocessing step. This further includes, Each preprocessing step included in the aforementioned composite preprocessing step is executed based on the parameters set for each preprocessing step. The computer implementation method according to claim 1.

9. The eddy current analysis mentioned above includes multiple independent conductor regions, The aforementioned multiple parts include a part that contains a discretized equation of the current continuity equation, The aforementioned composite preprocessing includes preprocessing for the part that includes a discretized equation of the current continuity equation, The pretreatment for the part including the discretized equation for current continuity comprises a plurality of pretreatments corresponding to the plurality of independent conductor regions, The plurality of preprocessing steps corresponding to the plurality of independent conductor regions are processed in parallel. The computer implementation method according to claim 1.

10. The aforementioned combined pretreatment includes a first pretreatment and a second pretreatment of a different type from the first pretreatment. The first preprocessing and the second preprocessing are performed in parallel. The computer implementation method according to claim 1.

11. The eddy current analysis mentioned above includes multiple independent conductor regions, The aforementioned parts are, Part 1 includes the discretized form of Ampere's Law, The second part includes the discretized equation for current continuity, Includes, The composite pretreatment includes a first pretreatment for the first part and a second pretreatment for the second part. The first preprocessing and the second preprocessing are performed in parallel. The second pretreatment comprises a plurality of pretreatments corresponding to the plurality of independent conductor regions, The plurality of preprocessing steps corresponding to the plurality of independent conductor regions are processed in parallel. The computer implementation method according to claim 1.

12. This system numerically analyzes eddy currents by solving the discretized governing equations using the finite element method with a pre-processed algorithm. The preprocessing in the aforementioned preprocessing algorithm is a composite preprocessing that performs different types of preprocessing on multiple parts of the governing equation. system.

13. This is a computer program that causes a computer to perform a process for numerically analyzing eddy currents by solving the discretized governing equations using the finite element method with a pre-processed algorithm. The preprocessing in the aforementioned preprocessing algorithm is a composite preprocessing that performs different types of preprocessing on multiple parts of the governing equation. Computer program.

Images