RECORDING MEDIUM RECORDING MAGNETIC MATERIAL SIMULATION PROGRAM, MAGNETIC MATERIAL SIMULATION METHOD AND INFORMATION PROCESSING APPARATUS

A recording medium recording a magnetic material simulation program includes: acquiring a shape model including element regions of a shape of a core, property information indicating physical properties of a magnetic material of the core, and coil current information indicating a time change in a current through a coil around the core; specifying a first current density of the coil at a first time, based on the coil current information; computing, using the property information and the first current density, first index values at the first time for positions in the shape model; computing, using the first index values, a charge density of each element region; specifying a second current density of the coil at a second time, based on the coil current information; and computing, using the property information, the second current density, and the charge density, second index values at the second time for the positions.

CROSS-REFERENCE TO RELATED APPLICATION

This application is based upon and claims the benefit of priority of the prior Japanese Patent Application No. 2017-030922, filed on Feb. 22, 2017, the entire contents of which are incorporated herein by reference.

FIELD

The embodiment discussed herein is related to a recording medium recording magnetic material simulation program, a magnetic material simulation method, and an information processing apparatus.

BACKGROUND

Along with improvements in the computational performance of computers, a variety of physical phenomena are being simulated using computers.

Related technology is disclosed in Japanese Laid-open Patent Publication No. 2006-351622.

SUMMARY

According to an aspect of the embodiments, a non-transitory computer-readable recording medium recording a magnetic material simulation program causing a computer to execute a process, the process includes: acquiring a shape model that includes a plurality of element regions generated by partitioning a shape of a core, property information indicating physical properties of a magnetic material used in the core, and coil current information indicating a time change in a current flowing through a coil wrapped around the core; specifying a first current density of the coil at a first time, based on the coil current information; computing, using the property information and the first current density, a plurality of first index values related to a magnetic field at the first time, in correspondence with a plurality of positions in the shape model; computing, using the plurality of first index values, a charge density of each of the plurality of element regions; specifying a second current density of the coil at a second time after the first time, based on the coil current information; and computing, using the property information, the second current density, and the charge density of each of the plurality of element regions, a plurality of second index values related to a magnetic field at the second time, in correspondence with the plurality of positions.

DESCRIPTION OF EMBODIMENT

For example, as one simulation of a physical phenomenon, for an electric component that includes a core using a magnetic material and a coil wrapped around the core, magnetic field analysis that analyzes the magnetic field produced when a current is made to flow through the coil is performed. For the magnetic field analysis simulation, the finite element method is used as the numerical analysis method.

With the finite element method, a model that partitions the shape of the analysis target into small regions called “elements” or a “mesh” is generated, and variables indicating unknown physical quantities are assigned to locations such as nodes, namely the centers of elements and vertices of elements, the edges of elements, and the like. A system of simultaneous equations (for example, a system of linear equations having a large-scale and sparse coefficient matrix) is generated from basic equations indicating the physical phenomenon and the multiple variables assigned in the model, and by solving the system of simultaneous equations, an approximate solution of the variables is obtained.

With a magnetic shield analysis method using the finite element method, an analysis routine including a first step, a second step, and a third step is executed. In the first step, a magnetic state is set for each of multiple sub-regions included in the shape of the analysis target. In the second step, the magnetic field states set in the first step and an integral equation are used to compute a boundary condition. In the third step, the boundary condition computed in the second step and the finite element method are used to update the magnetic field state of each sub-region. The above analysis routine is executed repeatedly while varying the time along the time axis. For example, the above three steps are executed again while accounting for the influence of eddy currents, and the magnetic field state of each sub-region is updated.

For example, with a core using a polycrystalline magnetic material, such as a manganese-zinc (MnZn) core, if a high-frequency current flows through a coil wrapped around the core, a dimensional resonance phenomenon may occur. The dimensional resonance phenomenon is a physical phenomenon in which the gaps existing between the multiple crystal grains included in a magnetic material behave like capacitors, thereby causing the permeability of the magnetic material to increase suddenly at specific frequencies. Since loss inside the core becomes large if the imaginary part of the permeability of the core increases, it is preferable to account for the dimensional resonance phenomenon in the magnetic field analysis of the core using a magnetic material. For example, in the magnetic field analysis of the finite element method, since the physical properties of the above magnetic material are not incorporated into the model or the basic equations, reproducing the dimensional resonance phenomenon by simulation may be difficult.

For example, a magnetic material simulation program and the like capable of improving the accuracy of magnetic field analysis of a magnetic material in which the dimensional resonance phenomenon occurs may be provided.

FIG. 1illustrates an exemplary magnetic material simulation apparatus. A magnetic material simulation apparatus10executes a magnetic field analysis of an electric component, which includes a core using a magnetic material and a coil wrapped around the core, as a computer simulation using the finite element method. The magnetic material simulation apparatus10may be client apparatus operated by a user, or a server apparatus.

The magnetic material simulation apparatus10includes a storage unit11and a processing unit12. The storage unit11is a volatile storage device such as random access memory (RAM), or a non-volatile storage device such as a hard disk drive (HDD) or flash memory. The processing unit12is a processor such as a central processing unit (CPU) or a digital signal processor (DSP), for example. The processing unit12may also include an electronic circuit for a specific purpose, such as an application-specific integrated circuit (ASIC) or a field-programmable gate array (FPGA). The processor executes programs stored in memory such as RAM. For example, a magnetic material simulation program is executed. Note that a set of multiple processors may be called a “multi-processor” or simply a “processor”.

The storage unit11stores a shape model13, property information14, and coil current information15used in the simulation. The shape model13, the property information14, and the coil current information15may be user-created, or generated by software such as modeling software. Additionally, the shape model13, the property information14, and the coil current information15may be created by the magnetic material simulation apparatus10, or acquired from another apparatus.

The shape model13is a finite element model that includes multiple element regions generated by partitioning the shape of the core. Examples of the core include donut-shaped and toroidal cores. The shape model13may additionally include element regions corresponding to the shape of the coil wrapped around the core, element regions corresponding to the space around the core and the coil, and the like. The multiple element regions are two-dimensional or three-dimensional regions of the same type, such as triangles, quadrilaterals, tetrahedrons, or hexahedrons, for example.

The property information14indicates the physical properties of the magnetic material used in the core. The property information14includes the conductivity and the permittivity of the magnetic material, for example. Also, the property information14includes, for example, size information indicating the relationship between the thickness of each of the crystal grains included in the magnetic material, and the thickness of the grain boundary layer existing between these multiple crystal grains. The coil current information15indicates change over time of the current flowing through the coil wrapped around the core. The change over time of the current may be expressed continuously by an amplitude or a frequency, or may be expressed discretely as a current density at each time.

The processing unit12specifies a current density16aof the current flowing through the coil at a time #1 in the simulation, based on the coil current information15. In the case in which the coil current information15is expressed by amplitude and frequency, the processing unit12obtains the current density16aat the time #1 by calculation. In the case in which the coil current information15is expressed discretely, the processing unit12selects the current density16aat the time #1 from among the current densities at multiple times.

The processing unit12uses the property information14and the current density16ato compute multiple index values related to the magnetic field at the time #1 in association with multiple positions on the shape model13. The multiple positions on the shape model13are positions specified based on the multiple element regions included in the shape model13, such as the centers of element regions, nodes enclosing the element regions, and edges enclosing the element regions. For example, the processing unit12computes an index value16bwith respect to an edge e1enclosing an element region, computes an index value16cwith respect to an edge e2, and computes an index value16dwith respect to an edge e3. The index values related to the magnetic field are numerical values indicating the physical phenomenon to be computed by the simulation, and are vector potentials or the like, for example. For example, from a designated calculation formula #1, the property information14, and the current density16a, the processing unit12generates a system of simultaneous equations including multiple variables assigned to multiple positions on the shape model13, and by solving the system of simultaneous equations, computes multiple index values as the solution of the multiple variables.

The processing unit12uses the multiple index values at the time #1, including the index values16b,16c, and16d, to compute the charge density of the charge accumulated in each of the multiple element regions included in the shape model13. For example, the processing unit12computes a charge density17of a certain element region using the index values16b,16c, and16dcorresponding to the edges e1, e2, and e3enclosing the element region. The charge density is computed in accordance with a designated calculation formula #2 different from the calculation formula #1 above, for example. At this time, the processing unit12may use the multiple index values at the time #1 to compute the current density of the current flowing through each of the multiple element regions included in the shape model13, and from the current density, compute the charge density of each of the multiple element regions.

The processing unit12specifies a current density18aof the current flowing through the coil at a time #2 after the time #1 in the simulation, based on the coil current information15. The processing unit12uses the property information14, the current density18a, and the charge density of each of the multiple element regions to compute multiple index values related to the magnetic field at the time #2 in association with multiple positions on the shape model13. For example, the processing unit12computes an index value18bwith respect to the edge e1, computes an index value18cwith respect to the edge e2, and computes an index value18dwith respect to the edge e3.

The index values18b,18c, and18dare of a similar type to the index values16b,16c, and16d, and are vector potentials or the like, for example. The calculation formula used at this point is basically similar to the calculation formula #1 above. However, the charge density17computed during the simulation at the time #1 is used in the calculation formula. For example, from the property information14, the current density18a, and the charge density of each of the multiple element regions, the processing unit12generates a system of simultaneous equations including the multiple variables above, and by solving the system of simultaneous equations, computes multiple index values as the solution of the multiple variables.

According to the magnetic material simulation apparatus10, in a magnetic material simulation using the finite element method, the charge density of the charge accumulated in each element region is computed from the simulation result at the time #1. The computed charge density is used in the simulation at the time #2 after the time #1 (for example, immediately after the time #1). With this arrangement, it is possible to reproduce the dimensional resonance phenomenon, in which the gaps existing between the multiple crystal grains included in the magnetic material behave like capacitors, thereby causing the permeability to increase suddenly at specific frequencies. Thus, the accuracy of the magnetic material simulation is improved.

For example, a simulation apparatus100executes, as a computer simulation, magnetic field analysis of an electric component including a MnZn toroidal core and a coil.

FIG. 2illustrates exemplary hardware of a simulation apparatus. The simulation apparatus100includes a CPU101, RAM102, an HDD103, an image signal processing unit104, an input signal processing unit105, a media reader106, and a communication interface107. The above units are coupled to a bus.

The CPU101is a processor including a computational circuit that executes the instructions of a program. The CPU101loads a program and at least partial data stored in the HDD103into the RAM102, and executes the program. Note that the CPU101may also be provided with multiple processor cores, the simulation apparatus100may also be provided with multiple processors, and the following processes may also be executed in parallel using multiple processors or processor cores. Also, a set of multiple processors may be called a “multi-processor” or a “processor”.

The RAM102is volatile semiconductor memory that temporarily stores programs executed by the CPU101and data used in computations by the CPU101. The simulation apparatus100may also be provided with a type of memory other than RAM, and may also be provided with multiple types of memory.

The HDD103is a non-volatile storage device that stores an operating system (OS), software programs such as application software, and data. The simulation apparatus100may be provided with another type of storage device, such as flash memory or a solid-state drive (SSD), and may also be provided with multiple non-volatile storage devices.

The image signal processing unit104, following instructions from the CPU101, outputs an image to a display31coupled to the simulation apparatus100. For the display31, a display of arbitrary type may be used, such as a cathode ray tube (CRT) display, a liquid crystal display (LCD), a plasma display, or an organic electro-luminescence (OEL) display.

The input signal processing unit105acquires an input signal from an input device32coupled to the simulation apparatus100, and outputs to the CPU101. For the input device32, a pointing device such as a mouse, a touch panel, a touchpad, or a trackball, a keyboard, a remote controller, buttons, switches, or the like may be used. Multiple types of input devices may be coupled to the simulation apparatus100.

The media reader106is a reading device that reads programs and data recorded onto a recording medium33. For the recording medium33, for example, a magnetic disk, an optical disc, a magneto-optical (MO) disc, semiconductor memory, or the like may be used. Magnetic disks include flexible disk (FD) and HDD. Optical discs include Compact Disc (CD) and Digital Versatile Disc (DVD).

For example, the media reader106copies programs and data read from the recording medium33to another recording medium, such as the RAM102or the HDD103. A read-out program is executed by the CPU101, for example. Note that the recording medium33may also be a portable recording medium, and may be used to distribute programs and data. The recording medium33and the HDD103are sometimes called computer-readable recording media.

The communication interface107is coupled to a network34, and is an interface that communicates with other apparatus via the network34. The communication interface107may be a wired communication interface coupled by a cable to a communication apparatus such as a switch, or a wireless communication interface coupled by a wireless link to a base station.

FIG. 3is a diagram illustrating an exemplary shape of a toroidal core and a coil. The electric component to be analyzed includes a toroidal core41and a coil42. In the toroidal core41, a polycrystalline magnetic material, namely a MnZn sintered magnet, is used. The sintered magnet is created by heating and compacting magnetic powder at high temperature. The toroidal core41has a donut shape with a central hole. The coil42is an electric wire wrapped around the toroidal core41. The coil42is wrapped spirally so as to pass through the hole of the toroidal core41. The coil42touches the surface of the toroidal core41, for example. Alternating current is supplied from the outside to the coil42. The current flowing through the coil42may be a high-frequency current.

The relative permeability of the toroidal core41(the ratio of the permeability of the toroidal core41with respect to the permeability of a vacuum) changes according to the frequency of the current flowing through the coil42. Specifically, if current of a specific, comparatively high frequency flows through the coil42, the dimensional resonance phenomenon, in which the relative permeability increases suddenly, occurs in the toroidal core41. The dimensional resonance phenomenon occurs due to the gaps (grain boundary layer) that exist between the multiple crystal grains included in the magnetic material behaving like capacitors.

FIG. 4illustrates an exemplary relationship between the frequency of the coil and the relative permeability of the toroidal core. The graph51illustrates the relationship between the frequency of current flowing through the coil42and the relative permeability of the toroidal core41. As illustrated in the graph51, each of the real part and the imaginary part of the relative permeability of the toroidal core41has a peak at a specific, comparatively high frequency. For example, the real part of the relative permeability takes a maximum value of approximately 1500 at approximately 8000 kHz. Also, for example, the imaginary part of the relative permeability takes a maximum value of approximately 1300 at approximately 9000 kHz.

If the imaginary part of the relative permeability of the toroidal core41becomes large, the loss per unit volume of the toroidal core41increases. For this reason, it is preferable for the designer of the electric component to ascertain the relationship between the frequency of the coil42and the loss of the toroidal core41. The simulation apparatus100performs magnetic field analysis reflecting such a dimensional resonance phenomenon.

FIG. 5illustrates an exemplary fine structure of a magnetic material. In the magnetic material used in the toroidal core41, there are formed multiple crystal grains, including a crystal grain61, and a grain boundary layer62existing between these multiple crystal grains. The crystal grain61has comparatively high conductivity, and comparatively low resistivity. The grain boundary layer62has comparatively low conductivity, that is, comparatively high resistivity. By enclosing the crystal grain61having high conductivity with the grain boundary layer62having low conductivity, eddy currents are suppressed. However, if a specific high-frequency current flows through the coil42, eddy currents occur readily in the crystal grain61.

Assume that inside the magnetic material, multiple crystal grains of uniform size are arranged regularly. When modeling the magnetic material, a single crystal grain and the grain boundary layer existing around the crystal grain are extracted as a single unit region. For the sake of simplicity, a two-dimensional unit region will be considered. As illustrated inFIG. 5, the crystal grain61and the grain boundary layer62enclosing three edges of the crystal grain61are extracted as a single unit region. The thickness of the grain boundary layer in a certain direction corresponds to the distance between adjacent crystal grains.

For such a unit region, potential differences V1and V2, current densities i1and i2, conductivities σ1and σ2, and distances L and D are defined. The potential difference V1is the difference in potential between one edge and the opposite edge of the crystal grain61. The potential difference V2is the difference in potential between one edge and the opposite edge of the grain boundary layer62(the edge of the crystal grain adjacent to one edge of the crystal grain61). The current density i1is the current density inside the crystal grain61. The current density i2is the current density inside the grain boundary layer62. The conductivity σ1is the conductivity of the crystal grain61. The conductivity σ2is the conductivity of the grain boundary layer62. The distance L is the distance between one edge and the opposite edge of the crystal grain61, or in other words, the thickness of the crystal grain61. The distance D is the distance between one edge and the opposite edge of the grain boundary layer62, or in other words, the thickness of the grain boundary layer62.

The electric properties of the unit region may be expressed by an equivalent electric circuit.FIG. 6is a diagram illustrating an exemplary electric circuit equivalent to the fine structure of the magnetic material. In the above unit region of the magnetic material, the grain boundary layer62may behave like a capacitor that accumulates charge. Accordingly, the electric properties of the above unit region is equivalent to the electric circuit ofFIG. 6. The electric circuit includes resistor circuits63and64, and a capacitor65. The resistor circuit64and the capacitor65are coupled in parallel to each other to the resistor circuit63. The resistor circuit63corresponds to the crystal grain61, while the resistor circuit64and the capacitor65correspond to the grain boundary layer62.

For example, the potential difference between one end and the other end of the resistor circuit63is V1. The potential difference between one end and the other end of the resistor circuit64is V2. Furthermore, for the electric circuit ofFIG. 6, resistances R1and R2, currents I2, I2, and IC, a charge Q, and a capacitance C are defined. The resistance R1is the resistance of the resistor circuit63, that is, the resistance of the crystal grain61(in units of ohms (Ω)). The resistance R2is the resistance of the resistor circuit64, that is, the resistance of the grain boundary layer62. The current I1is the current flowing through the resistor circuit63. The current I2is the current flowing through the resistor circuit64. The current ICis the current flowing through the capacitor65, and corresponds to I1-I2. The charge Q is the charge accumulated in the capacitor65, that is, the charge of the grain boundary layer62(in units of coulombs (C)). The capacitance C is the capacitance (electrostatic capacitance) of the capacitor65, that is, the capacitance of the grain boundary layer62(in units of farads (F)). Also, a charge density q is defined. The charge density q is the charge density of the capacitor65, that is, the charge density of the grain boundary layer62(in units of C/m2)).

Based on the above electric circuit model, the resistances R1and R2, the capacitance C, and the charge Q may be defined as in the Formulas (1) to (4). In Formula (3), εris the relative permittivity of the grain boundary layer62, while ε0is the permittivity of a vacuum.

The total potential V, which is the sum of the potential difference V1of the crystal grain61and the potential difference V2of the grain boundary layer62, is computed from the above Formulas (1) and (2), as in Formula (5). Provided that E is the average electric field of the unit region, since the distance D is sufficiently smaller than the distance L, the total potential V is approximated as in Formula (6). The average electric field E is computed from Formulas (5) and (6), as in Formula (7). In Formula (7), the dimensional ratio r is the ratio of the distance D of the grain boundary layer62with respect to the distance L of the crystal grain61.

Formula (8) is established as a formula of current conservation indicating the relationship among the currents I1, I2, and ICdescribed above, and from Formula (8) and Formula (4), Formula (9) is established as a formula of current density conservation.

Formula (10) is established from the relationship in which the potential difference across the resistor circuit64and the potential difference across the capacitor65are both V2. From Formula (10) and Formulas (2) to (4), Formula (11) is established between the current density i2and the charge density q. If Formula (11) is substituted into the current density i2of Formula (7), the average electric field E is computed as in Formula (12). If Formulas (12) is rearranged, the current density i1is computed as in Formula (13).

The charge density q of the charge accumulated in the grain boundary layer62may be computed by time-integrating the difference between the current density i1of the crystal grain61and the current density i2of the grain boundary layer62. Thus, utilizing Formula (11), the charge density q is computed as in Formula (14). In so doing, the relationship between the current density i1of the crystal grain61and the charge density q of the grain boundary layer62is defined.

In the case in which the capacitance effect of the grain boundary layer62is not considered, Formulas (15) and (16) may be used as the basic equations of magnetic field analysis. Formula (15) is a governing equation of a magnetic field, while Formula (16) is a formula of current density conservation. In Formulas (15) and (16), σ is the conductivity, A is the vector potential, and φ is the scalar potential. Also, in Formula (15), μ is the permeability, and J0is the current density of the coil42. Note that in the following, a variable that takes a vector value in a formula may be denoted with a vector symbol. The vector potential A of Formulas (15) and (16) as well as the current density J0of Formula (15) take vector values. In contrast, the scalar potential φ of Formulas (15) and (16) takes a scalar value. However, in the following, a certain variable may not be clearly indicated as taking a scalar value or a vector value.

The eddy current density i of the unit region is computed as in Formula (17), using the vector potential A and the scalar potential φ determined by Formulas (15) and (16). The eddy current density i takes a vector value. The eddy current density i is defined centered on the unit region, and is defined as the current density produced in accordance with the average electric field E. In other words, the relationship i=σE is established.

Next, the above basic equation is rearranged to account for the capacitance effect of the grain boundary layer62. Based on the fine structure ofFIG. 5, the crystal grain61is sufficiently larger than the grain boundary layer62, and the current density i1of the crystal grain61occupies the majority of the current density for the unit region overall. For this reason, the current density i1and the conductivity σ1of the crystal grain61appearing in Formula (13) may be treated as the current density and the conductivity representative of the unit region, and may be considered to be defined centered on the unit region. Accordingly, the eddy current density i of the unit region is approximated by the current density i1computed as in Formula (18). The current density i1depends on the vector potential A, the scalar potential φ, and the charge density q. Unlike Formula (17), a charge density q term has been added to Formula (18). The current i1and the charge density q are defined at the center position of the unit region.

Applying the current density i1and the conductivity σ1, the governing equation of the magnetic field of Formula (15) may be rearranged as in Formula (19). Unlike Formula (15), a charge density q term has been added to Formula (19). Also, the formula of current density conservation of Formula (16) may be rearranged as in Formula (20). Unlike Formula (16), a charge density q term has been added to Formula (20).

From Formula (14) described above, Formula (21) is established for the current density i1and the charge density q, each of which takes a vector value. If the time integral of Formula (21) is expressed as a numerical integration in the time direction, Formula (22) is established. In Formula (22), qnis the charge density at a time n inside the simulation, and qn+1is the charge density at a time n+1 (the next time after the time n) inside the simulation. In Formula (22), Δt is the minimum unit of the passage of time inside the simulation, and is the difference between the time n and the time n+1. In Formula (22), Δt has meaning as the time period over which to observe change in the charge density J0of the coil42. By rearranging Formula (22), Formula (23) is established. The charge density q at the time n+1 depends on the current density i1at the time n+1, and the charge density q at the immediately preceding time n.

The simulation apparatus100is able to use the above Formulas (18) to (20) and (23) as the basic equations of magnetic field analysis that account for the capacitance effect of the grain boundary layer62. For example, from the charge density J0at the time n+1 and the charge density q at the time n, the vector potential A at the time n+1 and the scalar potential φ at the time n+1 are computed in accordance with Formulas (19) and (20). Next, from the vector potential A at the time n+1, the scalar potential φ at the time n+1, and the charge density q at the time n, the current density i1at the time n+1 is computed in accordance with Formula (18). Subsequently, from the current density i1at the time n+1 and the charge density q at the time n, the charge density q at the time n+1 is computed in accordance with Formula (23). By repeatedly executing the above at multiple discrete times along the time axis, magnetic field analysis that accounts for the capacitance effect becomes possible.

FIG. 7illustrates an exemplary finite element model of a toroidal core and a coil. The shape of the electric component that includes the toroidal core41and the coil42is expressed by a model70. The model70is a set of elements, which are subdivided regions. The model70includes elements indicating the region of the toroidal core41, and elements indicating the region of the coil42. The model70may also include elements indicating air regions around the toroidal core41and the coil42.

Each element is a closed region enclosed by multiple nodes and multiple edges. The shape of each element is a triangle, a quadrilateral, a tetrahedron, a hexahedron, or the like, for example. A single triangular element is formed by three nodes and three edges. A single tetrahedral element is formed by four nodes and six edges. The same node may be shared by two or more elements. Also, the same edge may be shared by two or more elements.

FIG. 8illustrates an exemplary assignment of variables to a finite element model. Herein, for the sake of simplicity, a two-dimensional triangular element71will be considered. The element71is formed by nodes72,73, and74, and by edges75,76, and77. The edge75is a line segment joining the node72and the node73. The edge76is a line segment joining the node73and the node74. The edge77is a line segment joining the node74and the node72.

If given the model70to use in magnetic field analysis, the simulation apparatus100assigns an unknown vector potential A to each of the multiple edges included in the model70. Also, the simulation apparatus100assigns an unknown scalar potential φ to each of the multiple nodes included in the model70. The unknown vector potentials A and the unknown scalar potentials φ correspond to the variables of a system of simultaneous equations. Each of the assigned vector potentials A and scalar potentials φ takes a vector value and a scalar value, respectively.

Through simulation, an approximate solution to the vector potentials A and the scalar potentials φ is computed. The current density i1that takes a vector value and the charge density q that takes a vector value are assigned to the center position of each element. As described earlier, the current density i1and the charge density q depend on the vector potential A and the scalar potential φ. The current density i1and the charge density q of a certain element are computed from the scalar potentials φ assigned to the nodes forming the element, and the vector potentials A assigned to the edges forming the element.

For example, an unknown scalar potential φ1is assigned to the node72. An unknown scalar potential φ2is assigned to the node73. An unknown scalar potential φ3is assigned to the node74. Also, an unknown vector potential A1is assigned to the edge75. An unknown vector potential A2is assigned to the edge76. An unknown vector potential A3is assigned to the edge77. Also, the current density i1and the charge density q are assigned to the center of the element71. The current density i1and the charge density q depend on A1, A2, and A3, and also on φ1, φ2, and φ3.

A method of generating a system of simultaneous equations will be described. The vector potential of a certain element is computed as in Formula (24). In Formula (24), A on the left-hand side is a vector potential defined at the center of the element, and takes a vector value. Also, neis the number of edges forming a single element. In the case in which the element is a two-dimensional triangle, ne=3, and in the case in which the element is a three-dimensional tetrahedron, ne=6. Aiis the unknown vector potential assigned to the ith edge, and takes a scalar value. Neiis an interpolation function with respect to the ith edge, and has vector-like properties. The vector potential A of each element is computed from Ai, the number of which is ne, and from Nei, the number of which is ne.

The scalar potential of a certain element is computed as in Formula (25). In Formula (25), φ on the left-hand side is a scalar potential defined at the center of the element, and takes a scalar value. Also, nnis the number of nodes forming a single element. In the case in which the element is a two-dimensional triangle, nn=3, and in the case in which the element is a three-dimensional tetrahedron, nn=4. Also, φiis the unknown scalar potential assigned to the ith node, and takes a scalar value. Nniis an interpolation function with respect to the ith node, and has scalar-like properties. The scalar potential φ of each element is computed from φi, the number of which is nn, and from Nni, the number of which is N.

In Formula (19), which is one of the basic equations, Formula (24) is substituted in for the vector potential A, and Formula (25) is substituted in for the scalar potential φ. Subsequently, both sides of Formula (19) are multiplied by a weighting function w defined by Formula (26), and spatially integrated. In Formula (26), Neiis an interpolation function similar to Formula (24), while weiis a scalar value expressing the weight of an edge. By the above calculation, Formula (27) is computed. From Formula (27), there is generated a system of simultaneous equations including a number of equations corresponding to the number of vector potentials Ai, that is, the number of edges included in the model70.

The first term on the left side of Formula (27) may be expanded as in Formula (28) using matrices. Herein, {we} is a matrix (row vector) of the size 1×(number of edges), [Ne·Ne] is a matrix of the size (number of edges)×(number of edges), and {An+1-An} is a matrix (column vector) of the size (number of edges)×1. The second term on the left side of Formula (27) may be expanded as in Formula (29) using matrices. Herein, [Ne·∇Nn] is a matrix of the size (number of edges)×(number of nodes), and {φn+1} is a matrix of the size (number of nodes)×1. The third term on the left side of Formula (27) may be expanded as in Formula (30) using a matrix. Herein, [Ne·q] is a matrix of the size (number of edges)×1. The fourth term on the left side of Formula (27) may be expanded as in Formula (31) using matrices. Herein, [V×Ne·V×Ne] is a matrix of the size (number of edges)×(number of edges), and {An+1} is a matrix of the size (number of edges)×1. The right side of Formula (27) may be expanded as in Formula (32) using a matrix. Herein, [Ne·J0] is a matrix of the size (number of edges)×1.

If Formulas (27) to (32) are combined, Formula (33) is obtained. From Formula (33), there is generated a system of simultaneous equations including a number of equations corresponding to the number of edges included in the model70.

Similarly, in Formula (20), which is one of the basic equations, Formula (24) is substituted in for the vector potential A, and Formula (25) is substituted in for the scalar potential φ. Both sides of Formula (20) are multiplied by a weighting function N defined by Formula (34), and spatially integrated. In Formula (34), Nniis an interpolation function similar to Formula (25), while wniis a scalar value expressing the weight of a node.

The first term on the left side of the formula computed by the above calculation may be expanded as in Formula (35) using matrices. Herein, {wn} is a matrix of the size 1×(number of nodes), [∇Nn·Ne] is a matrix of the size (number of nodes)×(number of edges), and {An+1-An} is a matrix of the size (number of edges)×1. The second term on the left side of the computed formula may be expanded as in Formula (36) using matrices. Herein, [∇Nn·∇Nn] is a matrix of the size (number of nodes)×(number of nodes), and {φn+1} is a matrix of the size (number of nodes)×1. The third term on the left side of the computed formula may be expanded as in Formula (37) using a matrix. Herein, [∇Nn·q] is a matrix of the size (number of nodes)×1.

If Formulas (35) to (37) are combined, Formula (38) is obtained. From Formula (38), there is generated a system of simultaneous equations including a number of equations corresponding to the number of nodes included in the model70.

If the above Formulas (35) and (38) are united and combined, Formula (45) is obtained as the final system of simultaneous equations. Formula (45) expresses that the product of a coefficient matrix G and a solution vector X is equal to a right-hand side vector F. The coefficient matrix G is the union of the four partial coefficient matrices GAA, GAφ, GφA, and Gφφ. The partial coefficient matrix GAAis a matrix of the size (number of edges)×(number of edges) defined by Formula (39). The partial coefficient matrix GAφis a matrix of the size (number of edges)×(number of nodes) defined by Formula (40). The partial coefficient matrix GφAis a matrix of the size (number of nodes)×(number of edges) defined by Formula (41). The partial coefficient matrix Gφφis a matrix of the size (number of nodes)×(number of nodes) defined by Formula (42).

In the coefficient matrix G, GAAis disposed in the upper-left, GAφis disposed in the upper-right, GφAis disposed in the lower-left, and Gφφis disposed in the lower-right. The coefficient matrix G is a square matrix in which the length of one edge is (number of edges)+(number of nodes). The solution vector X is a column vector in which the unknown vector potentials A assigned to the edges of the model70and the unknown scalar potentials φ assigned to the nodes of the model70are arranged. The solution vector X corresponds to a variable vector in which the variables of the system of simultaneous equations are arranged. The solution vector X is a column vector in which the number of rows is (number of edges)+(number of nodes).

The right-hand side vector F is the union of two partial vectors FAand Fφ. The partial vector FAis a column vector of the size (number of edges)×1 defined by Formula (43). The partial vector FAis computed using the vector potential A and the charge density q of the preceding time. The partial vector Fφis a column vector of the size (number of nodes)×1 defined by Formula (44). The partial vector Fφis computed using the vector potential A and the charge density q of the preceding time. In the right-hand side vector F, FAis disposed on the top, and Fφis disposed on the bottom. The right-hand side vector F is a column vector in which the number of rows is (number of edges)+(number of nodes).

In this way, in accordance with the basic equations, namely Formulas (19) and (20), there is generated a system of simultaneous equations including equations corresponding to the sum of the number of edges and the number of nodes included in the model70, for example, the sum of the number of unknown vector potentials and the number of unknown scalar potentials. The solution to the system of simultaneous equations indicated by Formula (45) may be computed by an iterative method such as the conjugate gradient (CG) method or the incomplete Cholesky CG (ICCG) method.

Formula (18) indicating the current density i1may be transformed as in Formula (46) in a discrete time simulation. From the vector potential An+1at the time n+1, the scalar potential φn+1at the time n+1, the vector potential Anat the time n, and the charge density qnat the time n, the current density i1n+1of each element at the time n+1 is computed in accordance with Formula (46). From the current density i1n+1at the time n+1 and the charge density qnat the time n, the charge density qn+1of each element at the time n+1 is computed in accordance with Formula (23) described earlier.

FIG. 9illustrates exemplary functions of a simulation apparatus. The simulation apparatus100includes a model storage unit111, a parameter storage unit112, a result storage unit113, an equation generation unit121, a solution computation unit122, and a result display unit123. The model storage unit111, the parameter storage unit112, and the result storage unit113are implemented using a storage area reserved in the RAM102or the HDD103, for example. The equation generation unit121, the solution computation unit122, and the result display unit123are implemented using a program module executed by the CPU101, for example.

The model storage unit111stores model data indicating the multiple elements included in the model70. The model data may be generated based on user operations, or may be generated automatically by software such as modeling software. The model data may be generated by the simulation apparatus100, or may be generated by another apparatus.

The result storage unit113stores simulation results. The simulation results include vector potentials, scalar potentials, current densities, and charge densities. The equation generation unit121assigns variables to multiple positions in the model70, based on the model data stored in the model storage unit111. The equation generation unit121assigns vector potential variables to the edges of the elements, and assigns scalar potential variables to the nodes of the elements. Additionally, the equation generation unit121uses the parameters stored in the parameter storage unit112and the simulation results of the preceding step stored in the result storage unit113to generate a system of simultaneous equations that obtains a solution to the variables. The equation generation unit121repeats the generation of the system of simultaneous equations along the time axis, in accordance with the parameters.

The solution computation unit122uses an iterative method to solve the system of simultaneous equations generated by the equation generation unit121, and obtains an approximate solution to the variables. The solution computation unit122obtains an approximate solution to the vector potentials of the edges of the elements, and an approximate solution to the scalar potentials of the nodes of the elements. Also, the solution computation unit122uses the most recent vector potentials and scalar potentials as well as the simulation results of the preceding step stored in the result storage unit113to compute the current density and the charge density of each element. The solution computation unit122saves the vector potentials, the scalar potentials, the current densities, and the charge densities of the most recent step in the result storage unit113.

When the iterative processes by the equation generation unit121and the solution computation unit122end, the result display unit123causes the display31to display the simulation results stored in the result storage unit113. At this time, the result display unit123may also cause the display31to display numerical values included in the simulation results. In addition, the result display unit123may also cause the display31to display a visualization image converted from the simulation results. The result display unit123may also output the simulation results to an output apparatus other than the display31.

FIG. 10illustrates an exemplary node table and an element table. A node table114is stored in the model storage unit111. The node table114includes node number, X coordinate, and Y coordinate fields. The node number is an identification number of a node. The X coordinate and the Y coordinate are coordinates indicating the position of a node in an orthogonal coordinate system. In the case in which the model70is a three-dimensional mode, the node table114additionally includes a Z coordinate field.

A element table115is stored in the model storage unit111. The element table115includes element number, Node 1, Node 2, and Node 3 fields. The element number is an identification number of an element. Node 1, Node 2, and Node 3 indicate the node numbers of nodes forming an element. An edge that forms the element may be identified by a pair of two nodes forming the element. In the case in which the element is a tetrahedron, the element table115additionally includes a Node 4 field. For example, node numbers equal to the number of nodes forming an element are associated with a single element number.

FIG. 11illustrates an exemplary parameter table. The parameter table116is stored in the parameter storage unit112. The parameter table116associates parameter names and their values.

The parameters include a total number of steps and a time width Δt. The total number of steps is the number of discrete times in the simulation, and for example, expresses the number of repetitions of the process of generating a system of simultaneous equations and obtaining a solution. The time width Δt is the time difference in the simulation between a certain step and the next step, or in other words, expresses the minimum unit of the passage of time.

The parameters include a coil current and a coil current frequency. The coil current expresses the amplitude of the current flowing through the coil42. The coil current frequency expresses the frequency of the current flowing through the coil42. Assuming that the current flowing through the coil42varies sinusoidally, the current density J0of the coil42at each time may be computed from the coil current and the coil current frequency. The current density J0of each time may also be recorded in the parameter table116.

The parameters include physical property values for each of multiple materials. Material1is the air existing around the toroidal core41and the coil42. The parameters related to air include the permeability μ. Material2is the material of the coil42. The parameters related to the coil42include the permeability μ and the conductivity σ. Material3is the magnetic material used in the toroidal core41. The parameters related to the magnetic material include the permeability μ, the conductivity σ1, the relative permittivity εr, and the dimensional ratio r. Although the parameters of three materials are stated, the parameters of one or more materials are stated in the parameter table116, in accordance with the analysis target.

FIG. 12is a diagram illustrating an exemplary result table. The result table117is stored in the result storage unit113. The result table117includes step number, vector potential, scalar potential, current density, and charge density fields. The step number is a number that identifies each step (each time) of the simulation, and is an integer equal to or greater than 0, and less than or equal to the total number of steps.

The vector potential field lists the vector potentials of the multiple edges computed in a certain step. The vector potential of each edge in step 0 is initialized to 0. The scalar potential field lists the scalar potentials of the multiple nodes computed in a certain step. The scalar potentials in step 0 do not have to be computed. The current density field lists the current densities of the multiple elements computed in a certain step. The current densities in step 0 do not have to be computed. The charge density field lists the charge densities of the multiple elements computed in a certain step. However, the charge density of each element in step 0 is initialized to 0.

FIG. 13illustrates an exemplary procedure of a simulation.

(S10) The equation generation unit121reads out model data from the model storage unit111. The equation generation unit121reads out parameters from the parameter storage unit112.

(S11) The equation generation unit121assigns a vector potential variable (a variable indicating an unknown vector potential A) to each of the multiple edges indicated by the model data. The equation generation unit121assigns a scalar potential variable (a variable indicating an unknown scalar potential φ) to each of the multiple nodes indicated by the model data.

(S12) The equation generation unit121initializes the step number n to 0. The equation generation unit121specifies the total number of steps (an upper limit on the number of steps) indicated by the parameters.

(S13) The solution computation unit122initializes the value (vector potential A0) in step 0 of the vector potential variable assigned in step S11to 0. The solution computation unit122initializes the charge density q0in step 0 of each of the multiple elements indicated by the model data to 0. The solution computation unit122registers the vector potential A0and the charge density q0in the result table117stored in the result storage unit113.

(S14) The equation generation unit121computes the current density J0of the coil42in step n+1, based on the parameters read out in step S10.

(S15) The equation generation unit121, following Formula (45), generates the coefficient matrix G and the right-hand side vector F indicating the system of simultaneous equations corresponding to step n+1. In the generation of the coefficient matrix G, the parameters read out in step S10are used. In the generation of the right-hand side vector F, the current density J0computed in step S14and the parameters read out in step S10are used. Additionally, in the generation of the right-hand side vector F, the vector potential Anand the charge density qnof step n stored in the result table117are used.

(S16) The solution computation unit122computes a solution to the system of simultaneous equations generated in step S15with an iterative method such as the CG method or the ICCG method. With this arrangement, the vector potential An+1and the scalar potential φn+1of step n+1 are computed. The solution computation unit122registers the vector potential An+1and the scalar potential φn+1in the result table117.

(S17) The solution computation unit122, following Formula (46), computes the current density i1n+1of each element in step n+1. In the computation of the current density i1n+1, the read-out parameters, the vector potential An+1of step n+1, the vector potential Anof step n, the scalar potential φn+1of step n+1, and the charge density qnof step n are used. The solution computation unit122registers the current density i1n+1in the result table117.

(S18) The solution computation unit122, following Formula (23), computes the charge density qn+1of each element in step n+1. In the computation of the charge density qn+1, the read-out parameters, the current density i111+1of step n+1, and the charge density qnof step n are used. The solution computation unit122registers the charge density qn+1in the result table117.

(S19) The equation generation unit121determines whether the current step number n is less than the total number of steps. In the case in which the step number n is less than the total number of steps, the process proceeds to step S20, whereas in the case in which the step number n has reached the total number of steps, the simulation ends.

(S20) The equation generation unit121increments the step number n by 1, and proceeds to step S14.FIG. 14illustrates an exemplary simulation that does not account for a dimensional resonance phenomenon.

In the case of performing a simulation that does not account for the capacitance effect of the grain boundary layer62, in which the dimensional resonance phenomenon is not reproduced, a simulation result like the graph52is obtained. The simulation that does not account for the capacitance effect is a simulation using Formulas (15) to (17) described above as the basic equations. Similarly to the graph51described above, the graph52illustrates the relationship between the frequency of current flowing through the coil42and the relative permeability of the toroidal core41. Unlike the actual relative permeability of the toroidal core41illustrated by the graph51, the graph52does not express the increase in the real part and the increase in the imaginary part of the relative permeability, and does not reproduce the dimensional resonance phenomenon.

FIG. 15illustrates an exemplary simulation that accounts for the dimensional resonance phenomenon. In the case of performing a simulation that accounts for the capacitance effect of the grain boundary layer62, in which the dimensional resonance phenomenon is reproduced, a simulation result like the graph53is obtained. The simulation that accounts for the capacitance effect is a simulation using Formulas (18) to (20) and (23) described above as the basic equations. Similarly to the graphs51and52described above, the graph53illustrates the relationship between the frequency of current flowing through the coil42and the relative permeability of the toroidal core41. Similarly to the actual relative permeability of the illustrated by the graph51, the graph53expresses the increase in the real part and the increase in the imaginary part of the relative permeability, and reproduces the dimensional resonance phenomenon. In the graph53, the peak frequency at which the real part of the relative permeability reaches a maximum and the peak frequency at which the imaginary part reaches a maximum sufficiently approximate the actual peak frequencies illustrated by the graph51.

FIG. 16illustrates an exemplary magnetic flux density visualization image. As described above, the result display unit123is capable of converting result data registered in the result table117into a visualization image, and display the visualization image on the display31. An exemplary visualization image which may be generated by the result display unit123is a magnetic flux density visualization image54like inFIG. 16. The magnetic flux density visualization image54uses arrows to express the magnitude and direction of the magnetic flux density inside the toroidal core41. The magnetic flux density may be computed from the vector potentials A.

FIG. 17illustrates an exemplary current density visualization image. Another exemplary visualization image which may be generated by the result display unit123is a current density visualization image55like inFIG. 17. The current density visualization image55uses arrows to express the magnitude and direction of the eddy current density i inside the toroidal core41.

According to the simulation apparatus100illustrated inFIG. 2, in the magnetic material simulation using the finite element method, the capacitance effect of the grain boundary layer62is expressed as the charge density q of charge accumulated in each element, and the charge density q is added to the basic equations related to the magnetic field. The vector potential An+1and the scalar potential φn+1 at the time n+1 are computed from the vector potential Anand the charge density qnat the time n. The current density i1n(eddy current density) of each element at the time n+1 is computed from vector potential An+1and the scalar potential φn+1at the time n+1, and the vector potential Anand the charge density qnat the time n. The charge density qn+1of each element at the time n+1 is computed from the current density i1at the time n+1 and the charge density q at the time n, and is used in the calculation of the time n+2.

With this arrangement, the dimensional resonance phenomenon produced in a polycrystalline magnetic material may be reproduced, and the accuracy of the magnetic material simulation may be improved.