WO2019171660A1 - 磁場解析装置、解析方法、及びプログラム - Google Patents

磁場解析装置、解析方法、及びプログラム Download PDF

Info

Publication number
WO2019171660A1
WO2019171660A1 PCT/JP2018/041831 JP2018041831W WO2019171660A1 WO 2019171660 A1 WO2019171660 A1 WO 2019171660A1 JP 2018041831 W JP2018041831 W JP 2018041831W WO 2019171660 A1 WO2019171660 A1 WO 2019171660A1
Authority
WO
WIPO (PCT)
Prior art keywords
magnetic field
magnetization
volume element
volume
observation point
Prior art date
Application number
PCT/JP2018/041831
Other languages
English (en)
French (fr)
Inventor
修司 宮崎
Original Assignee
住友重機械工業株式会社
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 住友重機械工業株式会社 filed Critical 住友重機械工業株式会社
Priority to JP2020504776A priority Critical patent/JP7121111B2/ja
Publication of WO2019171660A1 publication Critical patent/WO2019171660A1/ja
Priority to US17/012,682 priority patent/US11703553B2/en

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/02—Measuring direction or magnitude of magnetic fields or magnetic flux
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/02—Measuring direction or magnitude of magnetic fields or magnetic flux
    • G01R33/10—Plotting field distribution ; Measuring field distribution
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N27/00—Investigating or analysing materials by the use of electric, electrochemical, or magnetic means
    • G01N27/72—Investigating or analysing materials by the use of electric, electrochemical, or magnetic means by investigating magnetic variables
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/12—Measuring magnetic properties of articles or specimens of solids or fluids
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16Z—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS, NOT OTHERWISE PROVIDED FOR
    • G16Z99/00—Subject matter not provided for in other main groups of this subclass
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/0023—Electronic aspects, e.g. circuits for stimulation, evaluation, control; Treating the measured signals; calibration

Definitions

  • the present invention relates to a magnetic field analysis device, an analysis method, and a program.
  • simulations incorporating magnetic field analysis are performed at the site of design and development of electric devices such as motors in accordance with the improvement of computer processing capacity.
  • simulation it is possible to evaluate to some extent without actually creating a prototype, so that the speed of design development can be improved.
  • a magnetic bead method has been proposed (for example, Patent Document 2).
  • the magnetic bead method an object is divided into a plurality of elements, each element is regarded as a spherical particle, and calculation is performed using an exact solution based on spherical symmetry.
  • the object is represented by an aggregate of a plurality of particles, a gap is generated between the particles, and therefore errors due to the gap and spherical particles may occur.
  • Patent Document 3 A method capable of shortening the calculation time of the magnetic field calculation has been proposed (for example, Patent Document 3).
  • a numerical solution of a magnetic field generated by a reference particle to which a magnetic moment is applied is calculated in advance for each of a plurality of lattice points set around the reference particle, and the vector amount of the magnetic moment applied to the reference particle is calculated. And the correspondence between the magnetic field and the numerical solution of the magnetic field are held for each lattice point.
  • a magnetic field generated by each particle to which a magnetic moment is applied is calculated using this correspondence. Since it is not necessary to calculate the area every time the magnetic field generated by each particle is calculated, the calculation time can be shortened.
  • JP 2009-37334 A Japanese Patent Laid-Open No. 2015-111401 JP 2017-194484 A
  • An object of the present invention is to provide a magnetic field analysis apparatus, an analysis method, and a program that can further reduce the calculation time.
  • a magnetization applying unit that divides the virtual space into a plurality of volume elements and applies magnetization to each of the volume elements; For each of a plurality of observation points in the virtual space, a magnetic field generated at the observation point for each volume element is calculated based on magnetization applied to the plurality of volume elements around the observation point, and a plurality of observation points are calculated.
  • a magnetic field calculation unit for obtaining a magnetic field generated at each of the plurality of observation points based on a calculation result for each volume element;
  • a magnetic field analysis device is provided in which the magnetization imparting unit makes a dimension in the first direction larger than a dimension in the other second direction for at least some of the volume elements.
  • the magnetic field applied to the volume element is generated for each of the volume elements to generate a magnetic field generated at a plurality of observation points in the virtual space, and based on the calculation results for the plurality of volume elements, Find the magnetic field generated in each
  • a magnetic field analysis method is provided in which at least a part of the volume elements has a dimension in the first direction in the virtual space larger than a dimension in the other second direction. .
  • a function of dividing a virtual space into a plurality of volume elements and imparting magnetization to each of the volume elements The magnetic field applied to the volume element is generated for each of the volume elements to generate a magnetic field generated at a plurality of observation points in the virtual space, and based on the calculation results for the plurality of volume elements, A function to determine the magnetic field generated by each,
  • the computer In dividing the virtual space into the volume elements, the computer is caused to realize a function of making at least some of the volume elements a size in the first direction in the virtual space larger than a size in the other second direction.
  • a program is provided.
  • the size of the volume element is smaller than the size of the volume element in the second direction, and the number of volume elements is smaller than the number of regular hexahedron voxels having the length in the second direction of the volume element as one side. . For this reason, compared with the case where it calculates for every voxel, the direction which calculates for every volume element becomes low. Further, if the direction in which the spatial change of the magnetic field is gentle is selected as the first direction, it is possible to suppress a decrease in calculation accuracy.
  • FIG. 1 is a block diagram of a magnetic field analysis apparatus according to this embodiment.
  • FIG. 2A is a perspective view of a state in which a part of the object is divided into a plurality of voxels
  • FIG. 2B is a perspective view of a part of the object after coarse graining.
  • FIG. 3A is a diagram showing a coordinate system defined in a virtual space
  • FIG. 3B is generated at the observation point by the magnetization of the volume element when the magnetization inside one volume element is assumed to be constant. it is a diagram for explaining a method of calculating the magnetic field H (r i).
  • FIG. 4 is a diagram illustrating a positional relationship between a volume element serving as a source point and a plurality of reference points around the volume element.
  • FIG. 4 is a diagram illustrating a positional relationship between a volume element serving as a source point and a plurality of reference points around the volume element.
  • FIG. 5 is a diagram showing the positional relationship between the voxel with magnetization and the observation point.
  • FIG. 6 is a diagram illustrating the positional relationship between the volume element serving as the source point and the observation point Po.
  • 7A and 7B are diagrams showing the positional relationship between the source point and the observation point.
  • 8A and 8B are a plan view and a side view, respectively, of a permanent magnet used in the magnet motor to be analyzed.
  • FIG. 8C is a result of magnetic field analysis performed by the method according to the embodiment and a magnetic field generated by the method according to the comparative example. It is a graph which shows the result of having analyzed.
  • FIG. 9 is a diagram showing a coordinate system used for explaining the method of multipole expansion.
  • FIG. 9 is a diagram showing a coordinate system used for explaining the method of multipole expansion.
  • FIG. 10A is a diagram for explaining a method of calculating a magnetic field H (r i ) generated at the observation point Po by performing dipole approximation for each voxel
  • FIG. 10B is an observation using a multipole expansion method. It is a figure explaining the method of calculating
  • FIG. 11 is a graph showing the error of the magnetic field Hc relative to the magnetic field Ho as a function of the normalized distance.
  • FIG. 12 is a graph showing the effect of speeding up the calculation by performing multipole expansion.
  • MD methods Molecular dynamics methods
  • Quantum quantum molecular dynamics methods
  • simulations based on Molecular Dynamics Method
  • RMD method renormalized molecular dynamics
  • the MD method and RMD method are kinetic methods (using statistical mechanics to calculate physical quantities)
  • the particle method is a method for discretizing differential equations describing continuums.
  • the MD method and the RMD method are also called particle methods here.
  • soot particle method can handle not only a static phenomenon but also a dynamic phenomenon such as a flow, it is attracting attention as a simulation method that replaces the finite element method that mainly analyzes a static phenomenon.
  • the soot particle method has a differential view of obtaining a particle system to be analyzed by discretizing a continuum with particles. For example, when dealing with fluids in the particle method, the Navier-Stokes equation is often discretized with particles.
  • another view of the particle method is an integral view of collecting many particles to form a continuum. This is, for example, the view that small iron grains are collected and hardened to form a large iron ball.
  • the magnetic fields generated at the observation points by the magnetic moments are added over the magnetic moments by the principle of superposition.
  • the present inventor has uniquely found the affinity between the method for obtaining a magnetic field from such a collection of magnetic moments and the particle method in the case of an integral view, and to give a magnetic moment to particles in the particle method. I came up with it. This makes it possible to extend the application range of the particle method to magnetic field analysis while maintaining the advantage of the particle method that is strong in convection and large deformation analysis.
  • the object is represented by an aggregate of spherical particles, and a high-speed and high-precision calculation result can be obtained by using an exact solution based on the spherical symmetry of each particle.
  • the object that is a continuum is represented by a collection of spherical particles, a gap is generated between the particles, and a calculation error due to the gap occurs.
  • the influence of the error becomes large.
  • FIG. 1 is a block diagram of a magnetic field analyzer 30 according to this embodiment.
  • the magnetic field analysis device 30 includes a processing device 20, a storage device 25, an input device 28, and an output device 29.
  • the processing device 20 includes an analysis information acquisition unit 21, a magnetization applying unit 22, a magnetic field calculation unit 23, and an output control unit 24.
  • the storage device 25 includes a data holding unit 26 and a correspondence relationship holding unit 27.
  • FIG. 1 can be realized in hardware by an element such as a central processing unit (CPU) or a mechanical device, and can be realized by a computer program or the like in software.
  • FIG. 1 shows functional blocks realized by cooperation of hardware and software. Therefore, these functional blocks can be realized in various modes by a combination of hardware and software.
  • the bag processing device 20 is connected to an input device 28 and an output device 29.
  • the input device 28 receives input of commands and data from a user related to processing executed by the processing device 20.
  • a keyboard or mouse that performs input by a user operation
  • a communication device that performs input via a network such as the Internet
  • a reading device that performs input from a recording medium such as a CD or DVD, or the like is used. Can do.
  • the analysis information acquisition unit 21 acquires magnetic field analysis information via the input device 28.
  • the magnetic field analysis information includes various information necessary for the magnetic field analysis. For example, the shape of the analysis object defined in the virtual space, the information for dividing the analysis object into voxels, the information for coarse-graining the voxels, and the spatial change of the magnetic field existing in the virtual space are relatively gradual Information for defining a specific direction is included.
  • the magnetization applying unit 22 divides the object in the virtual space into a plurality of voxels (regular hexahedral elements) based on the magnetic field analysis information.
  • Each voxel 50 is given a voxel identifier (voxel ID) for specifying one voxel 50.
  • FIG. 2A is a perspective view of a state in which a part of the object is divided into a plurality of voxels 50.
  • a plurality of voxels 50 are arranged in the x-axis direction, the y-axis direction, and the z-axis direction.
  • the plurality of voxels 50 have the same shape and dimensions, and the plurality of voxels 50 have a simple cubic unit cell structure.
  • the magnetic field includes an external magnetic field and a magnetic field generated by magnetizing an object.
  • the magnetization applying unit 22 applies magnetization to each of the voxels 50.
  • the magnetization may be given as an initial condition, or the magnetization generated in each voxel 50 may be calculated based on the magnetic field obtained by the magnetic field calculation unit 23 described later, and the magnetization obtained by the calculation may be given. Since the magnetization is a value per unit volume of the vector sum of all magnetic moments in the voxel 50, giving one value as the magnetization for each voxel 50 gives one value of the magnetic moment to the voxel 50. It is substantially the same as giving.
  • the magnetization applying unit 22 registers the applied magnetization in the data holding unit 26 in association with the voxel ID.
  • the magnetization imparting unit 22 collects a plurality of voxels 50 that are continuous in the direction in which the spatial change of the magnetic field is gentle (z-axis direction in the present embodiment) as one volume element 51.
  • the process of collecting a plurality of voxels 50 into one volume element 51 is called coarse graining.
  • the dimension of each of the coarse-grained volume elements 51 in the z-axis direction is larger than the dimensions in the x-axis direction and the y-axis direction.
  • FIG. 2B is a perspective view of a part of the target after coarse graining.
  • a plurality of voxels 50 (FIG. 2A) arranged in the z-axis direction are coarse-grained (collected) to form a volume element 51.
  • the number of voxels 50 coarse-grained into one volume element 51 may be set according to the degree of change in the magnetic field. For example, the number of voxels 50 coarse-grained into one volume element 51 may be increased as the change in the magnetic field is more gradual.
  • FIG. 2B shows an example in which the change in the z-axis direction of the magnetic field is gentle in the vicinity of the central portion in the vertical direction, and the change in the magnetic field gradually increases as it approaches the upper and lower ends.
  • the number of voxels 50 to be coarse-grained into one volume element 51 is relatively increased in the vicinity of the central portion in the vertical direction, and the coarser grained into one volume element 51 as it approaches the upper and lower ends.
  • the number of voxels 50 to be performed is relatively small.
  • the magnetization applying unit 22 applies magnetization to the volume elements 51 based on the magnetization applied to the original voxel 50. For example, a vector average of magnetization given to the original voxel 50 coarse-grained in one volume element 51 may be given to the volume element 51.
  • the magnetization applying unit 22 registers the magnetization applied to the volume element 51 in the data holding unit 26 in association with the volume element identifier (volume element ID).
  • the magnetic field calculation unit 23 calculates magnetic fields at a plurality of observation points in the virtual space based on the magnetization given to the voxel 50 or the magnetization given to the volume element 51, and holds the calculation results as data.
  • the voxel 50 or the volume element 51 having a magnetization that causes a magnetic field to be generated at the observation point is referred to as a source.
  • the observation point is arranged at the center position of each volume element 51, for example.
  • the magnetic field calculation unit 23 has three types of functions for calculating a magnetic field.
  • the first function (table reference method) is suitable for calculation when the distance from the source to the observation point is short, and uses the correspondence registered in the correspondence holding unit 27.
  • the second function (dipole approximation method) is suitable for calculation when the distance from the source to the observation point is long.
  • the third function (sequential integration method) is suitable for calculation when the observation point is close to the source. Details of the first to third functions will be described in detail later.
  • the output control unit 24 (FIG. 1) outputs the analysis result, for example, the calculation result of the magnetic field for each observation point registered in the data holding unit 26, to the output device 29.
  • FIG. 3A is a diagram showing a coordinate system defined in the virtual space. Represents the origin of the coordinate system in represents O, the position of the observation point Po in the virtual space at r i.
  • the position of the minute volume dV ′ is represented by r ′.
  • bold letters mean vectors.
  • the magnetic field H (r i ) at the observation point Po is obtained by integrating the magnetic field generated by the magnetization M (r ′) of the minute volume dV ′ over the entire space, and is expressed by the following equation.
  • ⁇ ′ means differentiation at the position r ′.
  • Equation (1) can be divided into terms of volume and area as follows.
  • dS ′ is a minute region on the surface of the region where the magnetic material is distributed
  • the vector ns is a unit normal vector of the surface S.
  • FIG. 3B illustrates a method of calculating the magnetic field H (r i ) generated at the observation point Po by the magnetization applied to the volume element 51 when the magnetization in one volume element 51 is constant.
  • the magnetization given to the j-th volume element 51 (position r j ) is represented by M j . It represents a small area of the surface of the volume element 51 in dS j, representing the unit normal vector of the micro area dS j at n s.
  • a magnetic field H (r i ) generated at the observation point Po at the position r i by transforming the expression (2) is expressed by the following equation.
  • r j represents the position of the minute region dS ′.
  • ⁇ regarding the parameter j means that the sum of the N volume elements 51 is taken.
  • the area for the surface S j means that the area is divided over the entire surface for each volume element 51.
  • FIG. 4 is a diagram showing the positional relationship between the volume elements 51S to be a source, a plurality of reference points G k around it.
  • FIG. 4 shows an example in which the volume element 51S as the source includes four continuous voxels 50.
  • it is generated in the magnetization M j of one volume element 51S serving as a source and a plurality of reference points G k (positions r i ).
  • the correspondence relationship with the magnetic field H k (r i ) to be calculated is calculated in advance by executing the area on the surface of the volume element 51S.
  • the reference point Gk is arranged corresponding to the positions of the vertices of the plurality of voxels 50, for example.
  • the magnetic field H k (r i ) generated at the reference point G k by the magnetization M j of the volume element 51S (position r j ) serving as the source can be expressed as follows from Equation (3).
  • the area means an area for the entire surface of the volume element 51S.
  • the magnetic field H k (r i ) represented by the equation (4) is an exact solution using the area.
  • the equation (4) is transformed and can be expressed as the following equation (5) using a tensor T k that defines the correspondence between the magnetic field H k (r i ) and the magnetization M j , any magnetization can be obtained.
  • the magnetic field H k (r i ) generated at each reference point G k can be easily calculated by the volume element 51 having M j .
  • n x, n y, n z is, x component of the unit normal vector n s, respectively, y component and z component.
  • x ij , y ij , and z ij are the x component, y component, and z component of the vector r ij starting from the minute region dS j of the volume element 51S serving as the source and ending at the reference point G k , respectively.
  • the volume element 51S serving as a source, the plurality of reference points G k running surface integral equation (6), are calculated in advance tensor T k for a plurality of reference points G k.
  • the tensor T k indicates the correspondence between the magnetization M j applied to the volume element 51S serving as the source and the magnetic field H k generated at each of the plurality of reference points G k by the magnetization M j .
  • This correspondence is prepared for each relative position between the volume element 51S serving as the source and the observation point.
  • the tensor T k indicating the correspondence is prepared for each size of the volume element 51 (FIG. 2B). For example, if the maximum number of voxels 50 to be coarse-grained is 10, tensor T k is prepared for each of the cases where the number of coarse-graining is 1 to 10.
  • the tensor T k indicating the correspondence relationship is stored in the correspondence relationship holding unit 27 (FIG. 1).
  • An example of a specific method for calculating the tensor T k in advance is described in Japanese Patent Application Laid-Open No. 2017-194484. Note that the method of executing the area of Equation (6) in advance is not limited to the example shown in Japanese Patent Application Laid-Open No. 2017-194484.
  • the magnetic field H (r i ) of the observation point Po (FIG. 3B) at the position r i is obtained by averaging the magnetic field H k of the reference point G k corresponding to the eight vertices of the voxel 50 including the observation point Po inside. Can be calculated.
  • the table reference method is used when the distance from the volume element 51S as the source to the observation point Po is equal to or less than a certain reference distance. Therefore, the tensor T k may be prepared for the reference point G k whose distance from the source volume element 51S is equal to or less than the reference distance.
  • the number of reference points G k for which the tensor T k should be prepared is smaller than the number of volume elements 51 in the virtual space. Therefore, the magnetic field is calculated using the tensor T k indicating the correspondence, compared to the method of calculating the magnetic field generated by all the volume elements 51 at the observation point Po (FIG. 3B) using the area for each time. The calculation load is lower when the table reference method is applied.
  • one volume element 51 after coarse graining becomes one source and one observation point. For this reason, if n voxels 50 are combined into one volume element 51, the number of sources and observation points are both 1 / n. Since the number of combinations of the source and the observation point is 1 / n 2 , the calculation time of the magnetic field decreases almost in proportion to n 2 .
  • FIG. 5 is a diagram showing a positional relationship between the voxel 50 provided with magnetization and the observation point Po.
  • the observation point Po is set at the center of each volume element 51, for example.
  • the magnetization from which the magnetic field is calculated the magnetization given to the voxel 50 before coarse graining is used.
  • the dipole approximation method is equivalent to obtaining the magnetic field by representing the magnetization M j given to the voxel 50 before coarse graining by the magnetic dipole arranged at the center of the voxel 50. If the volume of the voxel 50 is represented by ⁇ V j , the magnetic moment can be represented as ⁇ V j M j . In the dipole approximation method, it is not necessary to execute the area for each of the volume elements 51 (FIG. 2A) serving as the source points, so that the calculation load is reduced.
  • the magnetic field is calculated using the magnetization M j given to each of the plurality of voxels 50 constituting the volume element 51, but the plurality of voxels 50 are located at the center position of the volume element 51 after coarse graining.
  • the magnetic dipole assigned to may be subjected to multipole expansion, and the magnetic field at the observation point Po may be calculated based on the multipole.
  • one observation point Po is arranged corresponding to one volume element 51 after coarse graining, and one voxel 50 before coarse graining becomes one source. For this reason, when n voxels 50 are collectively coarsened into one volume element 51, the number of observation points Po becomes 1 / n, and the number of sources remains unchanged. As a result, when the number of combinations of the source and the observation point becomes 1 / n, the calculation time of the magnetic field is reduced substantially in proportion to n.
  • the table reference method is applied to the calculation of the magnetic field at the observation point Po inside a sphere with a certain radius around the source
  • the dipole approximation method is applied to the calculation of the magnetic field at the observation point Po outside the sphere. Good.
  • FIG. 6 is a diagram illustrating a positional relationship between the volume element 51S serving as a source and the observation point Po.
  • the observation point Po is set at the center of each volume element 51, for example.
  • the magnetization of the position r ′ in the volume element 51S as the source is represented by M (r ′).
  • the magnetic field H (r i ) generated at the position of the observation point Po by the magnetization M (r ′) imparted to the volume element 51S serving as the source is expressed by Equation (2).
  • Equation (2) Can be calculated based on.
  • the volume of Expression (2) is executed for the minute volume element dV ′ in the volume element 51S as the source, and the area is executed for the minute region dS ′ on the surface of the volume element 51S as the source. .
  • the right side of the equation (2) Since the term of the volume is 0, the magnetic field H (r i ) can be obtained by executing the area on the right side of the equation (2).
  • the sequential integration method without averaging the magnetic field of the position of the reference point G k using the table reference method to calculate the magnetic field of the observation point Po directly. For this reason, the calculation accuracy of the magnetic field can be improved as compared with the table reference method. In particular, the effect of applying the sequential integration method is remarkable when the source point is close to the observation point.
  • the magnetization M j within the volume element becomes a source 51S (FIG. 4) of magnetic field is calculated under the condition that it is fixed.
  • the sequential integration method can also be applied when the magnetization in the volume element 51 changes spatially.
  • one volume element 51 after coarse graining becomes one source and one observation point. For this reason, when n voxels 50 are collectively coarsened into one volume element 51, the number of sources and observation points are both reduced to 1 / n. As a result, since the number of combinations of the source and the observation point becomes 1 / n 2 , the calculation time of the magnetic field decreases in proportion to n 2 .
  • FIG. 7A is a diagram illustrating a positional relationship between the source Ps and the observation points Po1 and Po2.
  • a table reference method is applied to the calculation of the magnetic field generated at the observation point Po1 existing in the first region R1 inside the sphere with the radius a1 centered on the source Ps. Thereby, the calculation load can be reduced while ensuring the calculation accuracy.
  • a dipole approximation method is applied to the calculation of the magnetic field generated at the observation point Po2 existing in the second region R2 outside the first region R1. Thereby, it is possible to further reduce the calculation load while ensuring the calculation accuracy.
  • FIG. 7B is a diagram showing a positional relationship between the source Ps and the observation points Po2, Po3, Po4.
  • the first region R1 is divided into a third region R3 inside the sphere having a radius a2 (a2 ⁇ a1) centered on the source Ps and a fourth region R4 outside.
  • a table reference method is applied to the calculation of the magnetic field generated at the observation point Po4 in the fourth region R4, as in the example of FIG. 7A.
  • a sequential integration method is applied to the calculation of the magnetic field generated at the observation point Po3 in the third region R3. Thereby, calculation accuracy can be raised more.
  • FIG. 8A and 8B are a plan view and a side view, respectively, of a permanent magnet used for a magnet motor to be analyzed.
  • Ten permanent magnets 60 are arranged at equal intervals in the circumferential direction along the side surface of the cylinder whose central axis is the z-axis.
  • the radial direction is defined as an x-axis and a y-axis.
  • Each flat cross section of the permanent magnet 60 has an arc shape and is magnetized in the radial direction.
  • the magnetization directions of two permanent magnets 60 adjacent in the circumferential direction are opposite to each other.
  • the permanent magnet 60 is divided into a plurality of voxels 50 (FIG. 2A).
  • FIG. 8B shows the boundary in the z-axis direction of the volume element 51 after coarse graining, and does not show the boundary in the x-axis direction and the y-axis direction.
  • the change in the z-axis direction of the magnetic field formed by the permanent magnet 60 is more gradual than the changes in the x-axis direction and the y-axis direction. For this reason, coarse graining is performed in which a plurality of voxels 50 (FIG. 2A) continuous in the z-axis direction are collectively used as a volume element 51 (FIG. 2B).
  • the change in the magnetic field is relatively gentle near the center in the z direction, and the change in the magnetic field becomes steep near both ends in the z axis direction.
  • the number of voxels 50 put together in one volume element 51 is relatively increased near the center of the permanent magnet 60 in the z-axis direction, and is relatively decreased near both ends in the z-axis direction.
  • the number of voxels 50 was 234,952
  • the number of volume elements 51 after coarse graining was 26,064.
  • a magnetic field analysis was also performed on a comparative example that was not coarse-grained, that is, an example in which all of the volume elements 51 were each composed of one voxel 50.
  • the radius a2 of the spherical surface serving as the boundary between the third region R3 and the fourth region R4 in FIG. 7B is set to be four times the length of one side of the voxel 50, and the fourth region R4 and the second region R2
  • the radius a1 of the spherical surface serving as the boundary between and the length of one side of the voxel 50 is set to 16 times.
  • FIG. 8C is a graph showing the result of magnetic field analysis performed by the method according to the example and the result of magnetic field analysis performed by the method according to the comparative example.
  • the horizontal axis represents the angle in the rotation direction around the z axis in the unit of “radians”, and the vertical axis represents the strength of the magnetic field in arbitrary units.
  • the graph of FIG. 8C shows the strength of the magnetic field on the surface of the permanent magnet 60 in the middle in the z-axis direction.
  • the solid line in the graph indicates the magnetic field analysis result according to the example, and the circle symbol indicates the magnetic field analysis result according to the comparative example.
  • the magnetic field analysis results of the example and the comparative example are in good agreement. This means that sufficient calculation accuracy is obtained even when coarse-grained.
  • the calculation time when the magnetic field analysis was performed by the method according to the comparative example was 33.9 times the calculation time when the magnetic field analysis was performed by the method according to the example. It was confirmed that the calculation load can be reduced and the magnetic field analysis time can be shortened by coarse-graining.
  • FIG. 9 is a diagram showing a coordinate system used for explaining the method of multipole expansion.
  • One volume element 51 is obtained by coarse-graining a plurality of cubic voxels 50 arranged in a line in a direction in which the change of the magnetic field is relatively gentle.
  • the number of a plurality of voxels 50 to be expanded in multipole is represented by N.
  • the position of the center of the coarse-grained volumetric element 51 represented by r c. And that the position r c of multipole expansion point C.
  • An xyz orthogonal coordinate system is defined in which the multipole expansion point C is the origin, and the direction in which the magnetic field change is relatively gentle is the z-axis direction.
  • a plurality of voxels 50 having magnetic dipoles to be subjected to multipole expansion at the multipole expansion point C are located on the z-axis.
  • the position of the j-th voxel 50 is represented by r j -r c
  • the position of the observation point Po is expressed by r i -r c.
  • the j-th granted magnetized voxel 50 expressed by M j represents the volume of the j-th voxel 50 in [Delta] V j.
  • the angle formed by the z axis and the vector r i -r c from the multipole expansion point C to the observation point Po is represented by ⁇ i .
  • a distance r ij between the j-th voxel 50 and the observation point Po is expressed by the following formula (8).
  • r ic and r jc are the distance from the multipole expansion point C to the observation point Po, and the distance from the multipole expansion point C to the j-th voxel 50, respectively.
  • P L (cos ⁇ i ) is an L-order Legendre polynomial with cos ⁇ i as an argument. Equation (8) is described in J. D. JACKSON, CLASSICAL ELECTRODYNAMICS THIRD EDITION, p. 103.
  • a magnetic field H (r i ) generated at a distant observation point Po by the j-th voxel 50 having the magnetic moment m j can be approximated by the following equation.
  • M j and ⁇ V j are the magnetization applied to the j-th voxel 50 and the volume of the j-th voxel 50, respectively.
  • the magnetic field created by the volume element 51 in which the N voxels 50 are coarse-grained at the distant observation point Po is obtained by adding the equation (14) to the N voxels 50.
  • Expression (14) is added to N voxels 50, Expression (7) is obtained.
  • Equation (10) and Equation (12) are substituted into Equation (14) for rearrangement yields the following equation.
  • the expressions in the first, second, and third stages in the braces represent the x component, the y component, and the z component, respectively.
  • m jz is the z component of the magnetic moment m j .
  • the coefficient s 1 L shown in Expression (17), s 2 L and s 3 L can be calculated at the stage where the voxel 50 is coarse-grained and the volume element 51 is defined.
  • r ij does not appear as a position vector, and only ric appears. Therefore, the magnetic field volume element 51 after coarse-grained make the observation point Po is enabled to be calculated using the position r c of the multipole expansion point C, and the position r i of the observation point Po.
  • the multipole expansion point C it is possible to shorten the calculation time.
  • the length d j of a side of the N plurality of voxels 50 is storm into one volume element 51 are equal, is applied respectively to the N voxels 50 magnetization It is preferable to set the simulation conditions so as to satisfy the condition that M j is equal.
  • the values of the coefficients s 1 L , s 2 L , and s 3 L are calculated in advance up to the expansion order of the Legendre polynomial used when actually performing the simulation calculation, and stored in the storage device 25 (FIG. 1). Good.
  • FIG. 10A is a diagram for explaining a method of calculating the magnetic field H (r i ) generated at the observation point Po by performing dipole approximation for each voxel 50.
  • 100 voxels 50 are arranged side by side in the z-axis direction.
  • Magnetization M is given to each of the voxels 50.
  • the magnetic field generated by each voxel 50 at the observation point Po was obtained, and the magnetic field H (r i ) was obtained by adding the 100 voxels 50 together.
  • the magnetic field obtained by this calculation is expressed as Ho.
  • FIG. 10B is a diagram for explaining a method for calculating the magnetic field H (r i ) generated at the observation point Po by using the multipole expansion method.
  • the volume element 51 was defined by coarse-graining 100 voxels 50 (FIG. 10A), and the magnetic moment of the 100 voxels 50 was multipole expanded at the multipole expansion point C that is the center of the volume element 51.
  • the length of one side of each of the cubic voxels 50 is represented by d.
  • the coarse-grained scale becomes 100d.
  • the coefficients s 1 L , s 2 L , and s 3 L were calculated in advance using Equation (17), and the magnetic field H (r i ) at the observation point Po was obtained using Equation (16).
  • the magnetic field obtained by this calculation is expressed as Hc.
  • the calculation of the magnetic field Hc was performed for three cases in which the expansion order of the Legendre polynomial was 0th order, 2nd order, and 4th order.
  • FIG. 11 is a graph showing the error of the magnetic field Hc relative to the magnetic field Ho as a function of the normalized distance.
  • the standardized distance is obtained by standardizing the distance from the multipole expansion point C to the observation point Po with the coarse-grained scale 100d.
  • the error is defined as (Hc ⁇ Ho) / Ho.
  • the thin solid line, the broken line, and the thick solid line shown in the graph of FIG. 11 indicate errors when the expansion order of the Legendre polynomial is 0th order, 2nd order, and 4th order, respectively.
  • the error decreases as the degree of expansion of the Legendre polynomial increases. It can also be seen that the error decreases as the normalized distance increases. As shown in FIGS. 7A and 7B, when the sequential integration method or the table reference method is used in a range where the distance from the source Ps to the observation point Po is short, the dipole approximation method is applied in a range somewhat distant from the source Ps. The For this reason, it is considered that sufficiently high accuracy can be obtained even if the expansion order of the Legendre polynomial is lowered.
  • FIG. 12 is a graph showing the effect of speeding up the calculation by performing multipole expansion.
  • the computation time of the application of the method (FIG. 10A) to perform the dipole approximation for each voxel 50 at T 0 the computation time of the application of the multipole expansion technique in T C, fast effect It is defined as T 0 / T C. For example, if the calculated time T C becomes half the computation time T 0, speed effect is doubled.
  • the calculation time can be shortened while maintaining the desired accuracy by applying the multipole expansion method.
  • the expansion order of the Legendre polynomial is lowered, the calculation time is shortened, but the accuracy is lowered.
  • the expansion order of the Legendre polynomial may be determined according to the required accuracy and the allowable calculation time.
  • the change in the magnetic field in the z-axis direction illustrated in FIG. 2B is more gradual than the change in the magnetic field in the x-axis direction and the y-axis direction.
  • This example corresponds to a case where the magnetized object has a rod-like or cylindrical shape that is long in the z-axis direction.
  • the change in the magnetic field in the x-axis direction and the y-axis direction becomes slower than the change in the magnetic field in the z-axis direction.
  • the dimensions of the volume element 51 in the x-axis direction and the y-axis direction may be larger than the dimensions in the z-axis direction.
  • a spherical bead inscribed in the voxel 50 (FIG. 2A) is arranged, a magnetization (magnetic moment) is imparted to the bead, and a motion of a magnetic material composed of an aggregate of beads is simulated ( It can also be applied to the magnetic bead method.
  • the magnetic bead method is described in JP-A-2015-111401 and JP-A-2017-194484.
  • the volume of the voxel 50 is expressed as ⁇ V j / ⁇ j using the magnetic bead filling rate ⁇ j . Therefore, when the method according to the present embodiment is applied to the magnetic bead method, ⁇ V j / ⁇ j may be used instead of ⁇ V j in Equation (7).
  • each voxel 50 (FIG. 2A) may be provided with magnetization generated by an external magnetic field, and the magnetic field by this magnetization may be obtained.
  • the analysis target is divided into voxels 50 (FIG. 2A), and a plurality of continuous voxels 50 are coarse-grained to form the volume element 51.
  • the analysis target is directly passed through the voxel 50 without passing through the voxel 50. You may divide into 51.
  • one observation point is set in one volume element 51 (FIG. 2B), but a plurality of observation points may be set in one volume element 51, or an observation point for each voxel 50. May be set.
  • a plurality of observation points are set in one volume element 51 and magnetization is generated based on a magnetic field generated at the observation point, a distribution of magnetization is generated in one volume element 51.
  • the magnetic field may be calculated using the magnetization M j obtained by averaging the magnetization distribution in one volume element 51.
  • the magnetization M j for each voxel 50 (FIG. 5) serving as a source may be obtained from the magnetic field at the observation point.
  • the sequential integration method is applied, the magnetic field at the observation point may be calculated based on the non-uniform magnetization M (r ′) in the volume element 51 serving as the source.

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Condensed Matter Physics & Semiconductors (AREA)
  • Engineering & Computer Science (AREA)
  • Chemical & Material Sciences (AREA)
  • Theoretical Computer Science (AREA)
  • General Health & Medical Sciences (AREA)
  • Electrochemistry (AREA)
  • Analytical Chemistry (AREA)
  • Biochemistry (AREA)
  • Health & Medical Sciences (AREA)
  • Immunology (AREA)
  • Pathology (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Chemical Kinetics & Catalysis (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Measuring Magnetic Variables (AREA)

Abstract

【課題】計算時間をより短縮することが可能な磁場解析装置を提供する。 【解決手段】磁化付与部が、仮想空間を複数の体積要素に分割し、体積要素の各々に磁化を付与する。磁場演算部が、仮想空間の複数の観測点の各々について、観測点の周りの複数の体積要素に付与された磁化に基づいて、体積要素ごとに観測点に発生する磁場を演算し、複数の体積要素ごとの演算結果に基づいて、複数の観測点の各々に発生する磁場を求める。磁化付与部は、少なくとも一部の体積要素について、第１方向に関する寸法を、他の第２方向に関する寸法より大きくする。

Description

磁場解析装置、解析方法、及びプログラム
　本発明は、磁場解析装置、解析方法、及びプログラムに関する。
　近年、コンピュータの処理能力向上に伴い、モータ等の電気機器の設計開発の現場で、磁場解析を取り入れたシミュレーションが行われる。シミュレーションを行うと、実際にプロトタイプを作製しなくてもある程度の評価が可能となるので、設計開発のスピードが向上し得る。
　古典力学や量子力学等を基にコンピュータを用いて物質科学全般の現象を探るための方法として、分子動力学法を、マクロスケールの系を扱えるように発展させた繰り込み群分子動力学法に基づくシミュレーションが知られている（例えば、特許文献１）。分子動力学法のような粒子法は、静的な現象だけではなく流れ等の動的な現象も取り扱えるため、主として静的な現象を解析対象とする有限要素法等に代わるシミュレーション手法として注目されている。
　例えば、各粒子に磁気モーメントを付与し、各粒子の球対称性に基づく厳密解を利用して磁気的な物理量を演算対象とすることで、比較的精度の高いシミュレーション結果を高速に得ることができる磁気ビーズ法が提案されている（例えば、特許文献２）。磁気ビーズ法では、対象物を複数の要素に分割し、各要素を球形状の粒子とみなして球対称性に基づく厳密解を利用して演算を行う。このとき、対象物を複数の粒子の集合体で表すと、粒子間に隙間が生じるため、その隙間や球形状の粒子に起因する誤差が生じ得る。計算精度を高めるためには、対象物を立方体等のボクセルで隙間なく分割し、ボクセルの形状に応じた厳密解に基づく面積分を実行する必要がある。この面積分を高精度で行おうとすると計算時間が膨大になってしまう。
　磁場演算の計算時間の短縮化を図ることが可能な手法が提案されている（例えば、特許文献３）。この手法では、磁気モーメントが付与される基準粒子が発生させる磁場の数値解を、基準粒子の周りに設定される複数の格子点ごとに予め算出し、基準粒子に付与される磁気モーメントのベクトル量と、磁場の数値解との対応関係を格子点ごとに保持しておく。磁気モーメントが付与される各粒子が発生させる磁場を、この対応関係を用いて演算する。各粒子が発生させる磁場を計算する度に面積分の計算を行う必要がないため、計算時間を短縮することができる。
特開２００９－３７３３４号公報 特開２０１５－１１１４０１号公報 特開２０１７－１９４８８４号公報
　計算時間をより短縮する手法が望まれている。本発明の目的は、計算時間をより短縮することが可能な磁場解析装置、解析方法、及びプログラムを提供することである。
　本発明の一観点によると、
　仮想空間を複数の体積要素に分割し、前記体積要素の各々に磁化を付与する磁化付与部と、
　前記仮想空間の複数の観測点の各々について、前記観測点の周りの複数の前記体積要素に付与された磁化に基づいて、前記体積要素ごとに前記観測点に発生する磁場を演算し、複数の前記体積要素ごとの演算結果に基づいて、複数の前記観測点の各々に発生する磁場を求める磁場演算部と
を有し、
　前記磁化付与部は、少なくとも一部の前記体積要素について、第１方向に関する寸法を、他の第２方向に関する寸法より大きくする磁場解析装置が提供される。
　本発明の他の観点によると、
　仮想空間を複数の体積要素に分割し、前記体積要素の各々に磁化を付与し、
　前記体積要素に付与された磁化が前記仮想空間内の複数の観測点に発生させる磁場を前記体積要素ごとに計算し、複数の前記体積要素ごとの計算結果に基づいて、複数の前記観測点の各々に発生する磁場を求め、
　前記仮想空間を前記体積要素に分割するにあたり、少なくとも一部の前記体積要素について、前記仮想空間内の第１方向に関する寸法を、他の第２方向に関する寸法より大きくする磁場解析方法が提供される。
　本発明のさらに他の観点によると、
　仮想空間を複数の体積要素に分割し、前記体積要素の各々に磁化を付与する機能と、
　前記体積要素に付与された磁化が前記仮想空間内の複数の観測点に発生させる磁場を前記体積要素ごとに計算し、複数の前記体積要素ごとの計算結果に基づいて、複数の前記観測点の各々に発生する磁場を求める機能と、
　前記仮想空間を前記体積要素に分割するにあたり、少なくとも一部の前記体積要素について、前記仮想空間内の第１方向に関する寸法を、他の第２方向に関する寸法より大きくする機能と
をコンピュータに実現させるためのプログラムが提供される。
　体積要素の第１方向に関する寸法を第２方向に関する寸法より大きくすること、体積要素の第２方向の寸法を一辺の長さとする正六面体のボクセルの個数より、体積要素の個数の方が少なくなる。このため、ボクセルごとに演算を行う場合に比べて、体積要素ごとに演算を行う方が、演算負荷が低くなる。また、磁場の空間的な変化が緩やかな方向を第１方向に選ぶと、計算精度の低下を抑制することができる。
図１は、本実施例による磁場解析装置のブロック図である。 図２Ａは、対象物の一部分を複数のボクセルに分割した状態の斜視図であり、図２Ｂは、粗視化後の対象物の一部分の斜視図である。 図３Ａは、仮想空間内に定義される座標系を示す図であり、図３Ｂは、１つの体積要素の内部の磁化が一定であるとしたときに、体積要素の磁化によって観測点に発生する磁場Ｈ（ｒｉ）を計算する方法を説明する図である。 図４は、ソース点となる体積要素と、その周りの複数の基準点との位置関係を示す図である。 図５は、磁化が付与されたボクセルと観測点との位置関係を示す図である。 図６は、ソース点となる体積要素と観測点Ｐｏとの位置関係を示す図である。 図７Ａ及び図７Ｂは、ソース点と観測点との位置関係を示す図である。 図８Ａ及び図８Ｂは、それぞれ解析対象の磁石電動機に用いられる永久磁石の平面図及び側面図であり、図８Ｃは、実施例による方法で磁場解析を行った結果、及び比較例による方法で磁場解析を行った結果を示すグラフである。 図９は、多重極展開の方法を説明するために用いる座標系を示す図である。 図１０Ａは、個々のボクセルについて双極子近似を行い、観測点Ｐｏに発生する磁場Ｈ（ｒｉ）を計算により求める方法を説明する図であり、図１０Ｂは、多重極展開手法を用い、観測点Ｐｏに発生する磁場Ｈ（ｒｉ）を計算により求める方法を説明する図である。 図１１は、磁場Ｈｏに対する磁場Ｈｃの誤差を、規格化距離の関数として示すグラフである。 図１２は、多重極展開を行うことによる計算の高速化効果を示すグラフである。
  古典力学や量子力学等を基に計算機を用いて物質科学全般の現象を探るための方法として、分子動力学法（Molecular Dynamics Method、以下ＭＤ法と称す。）や、量子分子動力学法（Quantum Molecular Dynamics Method）や、マクロスケールの系を扱えるようにＭＤ法を発展させた繰り込み群分子動力学法（Renormalized Molecular Dynamics、以下ＲＭＤ法と称す。）に基づくシミュレーションが知られている。正確にはＭＤ法やＲＭＤ法は運動論的手法（物理量の算出には統計力学を使う）であり、粒子法は、連続体を記述する微分方程式を離散化する手法であり、両者は別ものであるが、ここではＭＤ法やＲＭＤ法も粒子法と呼ぶ。
  粒子法は静的な現象だけでなく流れなどの動的な現象をも取り扱えるので、主に静的な現象を解析対象とする有限要素法などに代わるシミュレーション手法として注目されている。
  粒子法には、連続体を粒子で離散化することにより解析対象の粒子系を得るという微分的な見方がある。例えば、粒子法において流体を扱う場合、ナビエストークス（Navier-Stokes）方程式を粒子で離散化することが多い。
  一方、粒子法の別の見方として、多くの粒子を集めて連続体を形成するという積分的な見方もある。これは例えば、小さな鉄の粒を集めて固めて大きな鉄球を形成するという見方である。
  一般に、多くの磁気モーメントが存在する空間内のある点（観測点）の磁場を求める場合、重ね合わせの原理により、各磁気モーメントが観測点に発生させる磁場を磁気モーメントに亘って足し合わせる。本発明者は、このような磁気モーメントの集まりから磁場を求める手法と、積分的な見方をした場合の粒子法との親和性を独自に見い出し、粒子法における粒子に磁気モーメントを付与することに想到した。これにより、対流や大変形の解析に強いという粒子法の利点を維持しつつ、粒子法の適用範囲を磁場解析にまで広げることが可能となる。
　対象物を球形の粒子の集合体で表し、各粒子の球対称性に基づく厳密解を利用することで、高速かつ高精度の演算結果を得ることができる。しかしながら、連続体である対象物を球形の粒子の集合体で表すと粒子間に隙間が生じるため、その隙間に起因する計算誤差が生じてしまう。特に、粒子間の距離が近い場合の粒子近傍の磁場を計算する場合に誤差の影響が大きくなる。計算精度を高めるためには、対象物を立方体やボロノイ多面体といったボクセルに分割し、これらのボクセルにより隙間なく対象物を分割することが好ましい。この場合、磁場を求めるために、ボクセルごとに面積分を実行する必要がある。この面積分を実行するには、計算時間が膨大となる。
　以下に説明する実施例による磁場解析装置では、計算時間の短縮化を図ることが可能である。
　図１は、本実施例による磁場解析装置３０のブロック図である。本実施例による磁場解析装置３０は、処理装置２０、記憶装置２５、入力装置２８、及び出力装置２９を含む。処理装置２０は、解析情報取得部２１、磁化付与部２２、磁場演算部２３、及び出力制御部２４を含む。記憶装置２５は、データ保持部２６及び対応関係保持部２７を含む。
　図１に示す各ブロックは、ハードウェア的には、コンピュータの中央処理ユニット（ＣＰＵ）をはじめとする素子や機械装置で実現することができ、ソフトウェア的にはコンピュータプログラム等によって実現することができる。図１では、ハードウェア及びソフトウェアの連携によって実現される機能ブロックが示されている。従って、これらの機能ブロックは、ハードウェア及びソフトウェアの組み合わせによって、種々の態様で実現することが可能である。
  処理装置２０は入力装置２８及び出力装置２９と接続される。入力装置２８は、処理装置２０で実行される処理に関係するユーザからのコマンド及びデータの入力を受ける。入力装置２８として、例えばユーザが操作を行うことにより入力を行うキーボードやマウス、インターネット等のネットワークを介して入力を行う通信装置、ＣＤ、ＤＶＤ等の記録媒体から入力を行う読取装置等を用いることができる。
　解析情報取得部２１は、入力装置２８を介して磁場解析情報を取得する。磁場解析情報には、磁場解析に必要な種々の情報が含まれる。例えば、仮想空間内に定義される解析対象物の形状、解析対象物をボクセルに分割する情報、ボクセルを粗視化する情報、仮想空間内に存在する磁場の空間的な変化が相対的に緩やかな方向を定義する情報等が含まれる。
　磁化付与部２２は、磁場解析情報に基づき仮想空間内の対象物を複数のボクセル（正六面体要素）に分割する。ボクセル５０の各々に、１つのボクセル５０を特定するためのボクセル識別子（ボクセルＩＤ）が付与される。
　図２Ａは、対象物の一部分を複数のボクセル５０に分割した状態の斜視図である。複数のボクセル５０が、ｘ軸方向、ｙ軸方向、及びｚ軸方向に並んでいる。複数のボクセル５０の形状及び寸法は同一であり、複数のボクセル５０は、単純立方単位胞構造を持つ。本実施例において、ｚ軸方向に関する磁場の変化が相対的に緩やかであり、ｘ軸方向及びｙ軸方向に関する磁場の変化が相対的に急峻であるある例について説明する。ここで、磁場には、外部磁場、及び対象物が磁化されることによって発生する磁場が含まれる。
　さらに、磁化付与部２２（図１）は、ボクセル５０の各々に磁化を付与する。磁化は初期条件として付与してもよいし、後述する磁場演算部２３により求められた磁場に基づいて各ボクセル５０に発生する磁化を計算し、計算によって求められた磁化を付与してもよい。なお、磁化は、ボクセル５０内の全磁気モーメントのベクトル和の単位体積当たりの値であるから、ボクセル５０ごとに磁化として１つの値を付与することは、ボクセル５０に磁気モーメントの１つの値を付与することと実質的に同一である。磁化付与部２２は、付与された磁化を、ボクセルＩＤと関連付けてデータ保持部２６に登録する。
　さらに、磁化付与部２２は、磁場の空間的な変化が緩やかな方向（本実施例においてｚ軸方向）に連続する複数のボクセル５０をまとめて１つの体積要素５１とする。複数のボクセル５０をまとめて１つの体積要素５１とする処理を粗視化という。粗視化された体積要素５１の各々のｚ軸方向の寸法は、ｘ軸方向及びｙ軸方向の寸法より大きくなる。
　図２Ｂは、粗視化後の対象物の一部分の斜視図である。ｚ軸方向に並ぶ複数のボクセル５０（図２Ａ）が粗視化されて（まとめられて）体積要素５１が形成されている。１つの体積要素５１に粗視化するボクセル５０の個数は、磁場の変化の度合いに応じて設定するとよい。例えば、磁場の変化が緩やかであるほど、１つの体積要素５１に粗視化するボクセル５０の個数を多くするとよい。図２Ｂでは、上下方向の中央部の近傍で磁場のｚ軸方向の変化が緩やかであり、上下端に近づくに従って、磁場の変化が徐々に大きくなる例を示している。このため、上下方向の中央部の近傍において、１個の体積要素５１に粗視化するボクセル５０の個数を相対的に多くし、上下端に近づくほど、１個の体積要素５１に粗視化するボクセル５０の個数を相対的に少なくしている。
　磁化付与部２２（図１）は、元のボクセル５０に付与されている磁化に基づいて体積要素５１に、それぞれ磁化を付与する。例えば、１つの体積要素５１に粗視化された元のボクセル５０に付与されている磁化のベクトル平均を、体積要素５１に付与するとよい。磁化付与部２２は、体積要素５１に付与された磁化を、体積要素識別子（体積要素ＩＤ）と関連付けてデータ保持部２６に登録する。
　磁場演算部２３（図１）は、ボクセル５０に付与された磁化、または体積要素５１に付与された磁化に基づいて、仮想空間内の複数の観測点における磁場を計算し、計算結果をデータ保持部２６に登録する。観測点に磁場を生じさせる元となる磁化を持つボクセル５０または体積要素５１をソースということとする。観測点は、例えば体積要素５１の各々の中心位置に配置する。磁場演算部２３は、磁場を計算する３種類の機能を備えている。第１機能（テーブル参照方式）は、ソースから観測点までの距離が短いときの計算に適しており、対応関係保持部２７に登録されている対応関係を利用する。第２機能（双極子近似方式）は、ソースから観測点までの距離が長いときの計算に適している。第３機能（逐次積分方式）は、観測点がソースに近接しているときの計算に適している。第１機能～第３機能の詳細については、後に詳しく説明する。
　出力制御部２４（図１）は、解析結果、例えばデータ保持部２６に登録されている観測点ごとの磁場の計算結果を、出力装置２９に出力する。
　［第１機能（テーブル参照方式）］
　次に、図３Ａ～図４を参照して、磁場演算部２３の第１機能（テーブル参照方式）について説明する。
　図３Ａは、仮想空間内に定義される座標系を示す図である。座標系の原点をＯで表し、仮想空間内の観測点Ｐｏの位置をｒｉで表す。微小体積ｄＶ’の位置をｒ’で表す。ここで、太文字はベクトルを意味する。観測点Ｐｏの磁場Ｈ（ｒｉ）は、微小体積ｄＶ’の磁化Ｍ（ｒ’）によって発生する磁場を全空間に亘って積分することにより得られ、以下の式で表される。
Figure JPOXMLDOC01-appb-M000001
 　ここで、∇’は、位置ｒ’における微分を意味する。
　磁性材料が分布する領域をＶで表し、表面をＳで表すと、式（１）は、以下のように体積分と面積分の項に分けることができる。
Figure JPOXMLDOC01-appb-M000002
 　ここで、ｄＳ’は、磁性材料が分布する領域の表面の微小領域であり、ベクトルｎｓは、表面Ｓの単位法線ベクトルである。
　図３Ｂは、１つの体積要素５１の内部の磁化が一定であるとしたときに、体積要素５１に付与された磁化によって観測点Ｐｏに発生する磁場Ｈ（ｒｉ）を計算する方法を説明する図である。ｊ番目の体積要素５１（位置ｒｊ）に付与された磁化をＭｊで表す。体積要素５１の表面の微小領域をｄＳｊで表し、微小領域ｄＳｊの単位法線ベクトルをｎｓで表す。
　式（２）を変形して、位置ｒｉの観測点Ｐｏに発生する磁場Ｈ（ｒｉ）は、以下の式で表される。
Figure JPOXMLDOC01-appb-M000003
 　ここで、ｒｊは、微小領域ｄＳ’の位置を表す。パラメータｊに関するΣは、Ｎ個の体積要素５１についての和をとることを意味する。表面Ｓｊに関する面積分は、体積要素５１ごとに表面全域で面積分を行うことを意味する。
　図４は、ソースとなる体積要素５１Ｓと、その周りの複数の基準点Ｇｋとの位置関係を示す図である。図４では、ソースとなる体積要素５１Ｓが、連続する４個のボクセル５０を含む例を示している。本実施例においては、式（３）の面積分の演算処理を高速化させるために、ソースとなる１つの体積要素５１Ｓの磁化Ｍｊと、複数の基準点Ｇｋ（位置ｒｉ）に発生する磁場Ｈｋ（ｒｉ）との対応関係を、体積要素５１Ｓの表面における面積分を実行することにより予め算出しておく。基準点Ｇｋは、例えば複数のボクセル５０の頂点の位置に対応して配置される。
　ソースとなる体積要素５１Ｓ（位置ｒｊ）の磁化Ｍｊによって基準点Ｇｋに発生する磁場Ｈｋ（ｒｉ）は、式（３）から、以下のように表すことができる。
Figure JPOXMLDOC01-appb-M000004
 　ここで、面積分は、体積要素５１Ｓの表面の全域についての面積分を意味する。
　式（４）で表される磁場Ｈｋ（ｒｉ）は、面積分を利用した厳密解である。ここで、式（４）を変形し、磁場Ｈｋ（ｒｉ）と磁化Ｍｊとの対応関係を定めるテンソルＴｋを用いて下記式（５）のように表すことができれば、任意の磁化Ｍｊを持つ体積要素５１によって各基準点Ｇｋに発生する磁場Ｈｋ（ｒｉ）を容易に計算できるようになる。
Figure JPOXMLDOC01-appb-M000005
　テンソルＴｋの各要素は、以下のように表すことができる。
Figure JPOXMLDOC01-appb-M000006
 　ここで、ｎｘ、ｎｙ、ｎｚは、それぞれ単位法線ベクトルｎｓのｘ成分、ｙ成分、ｚ成分である。ｘｉｊ、ｙｉｊ、ｚｉｊは、それぞれソースとなる体積要素５１Ｓの微小領域ｄＳｊを始点とし基準点Ｇｋを終点とするベクトルｒｉｊのｘ成分、ｙ成分、ｚ成分である。
　本実施例では、ソースとなる体積要素５１Ｓに対して、複数の基準点Ｇｋについて式（６）の面積分を実行し、テンソルＴｋを複数の基準点Ｇｋについて予め算出している。テンソルＴｋは、ソースとなる体積要素５１Ｓに付与された磁化Ｍｊと、磁化Ｍｊによって複数の基準点Ｇｋの各々に発生する磁場Ｈｋとの対応関係を示している。この対応関係は、ソースとなる体積要素５１Ｓと観測点との相対位置ごとに準備されている。さらに、対応関係を示すテンソルＴｋは、体積要素５１（図２Ｂ）の大きさごとに準備されている。例えば、粗視化するボクセル５０の最大個数を１０個とすると、粗視化の個数が１個の場合から１０個の場合までのそれぞれについてテンソルＴｋが準備される。
　対応関係を示すテンソルＴｋは、対応関係保持部２７（図１）に記憶されている。テンソルＴｋを予め計算する具体的な方法の一例が、特開２０１７－１９４８８４号公報に説明されている。なお、式（６）の面積分を予め実行しておく方法は、特開２０１７－１９４８８４号公報に示された例に限定されない。
　位置ｒｉの観測点Ｐｏ（図３Ｂ）の磁場Ｈ（ｒｉ）は、観測点Ｐｏを内部に含むボクセル５０の８個の頂点に相当する基準点Ｇｋの磁場Ｈｋを平均化することにより算出することができる。
　テーブル参照方式は、ソースとなる体積要素５１Ｓから観測点Ｐｏまでの距離が、ある基準距離以下のときに利用される。従って、ソースとなる体積要素５１Ｓからの距離がこの基準距離以下の基準点ＧｋについてテンソルＴｋを準備しておけばよい。テンソルＴｋを準備しておくべき基準点Ｇｋの個数は、仮想空間内の体積要素５１の個数より少ない。このため、全ての体積要素５１が観測点Ｐｏ（図３Ｂ）に発生させる磁場を、その都度面積分を用いて計算する方式に比べて、対応関係を示すテンソルＴｋを利用して磁場を計算するテーブル参照方式を適用する方が、演算負荷が低くなる。
　テーブル参照方式では、粗視化後の１つの体積要素５１が１つのソースとなり、かつ１つの観測点となる。このため、ｎ個のボクセル５０をまとめて１つの体積要素５１とすると、ソース及び観測点の個数が共に１／ｎになる。ソースと観測点との組み合わせの個数は１／ｎ２になるため、磁場の演算時間がほぼｎ２に比例して減少する。
　［第２機能（双極子近似方式）］
　次に、図５を参照して、磁場演算部２３（図１）の第２機能（双極子近似方式）について説明する。
　図５は、磁化が付与されたボクセル５０と観測点Ｐｏとの位置関係を示す図である。観測点Ｐｏは、例えば体積要素５１の各々の中心に設定される。磁場を計算する元となる磁化として、粗視化前のボクセル５０に付与された磁化を用いる。
　式（１）を部分積分して近似を行うと、以下の式（７）が得られる。
Figure JPOXMLDOC01-appb-M000007
　双極子近似方式は、粗視化前のボクセル５０に付与された磁化Ｍｊを、ボクセル５０の中心に配置された磁気双極子で代表させて磁場を求めることと等価である。ボクセル５０の体積をΔＶｊで表すと、磁気モーメントはΔＶｊＭｊと表すことができる。双極子近似方式においては、ソース点となる体積要素５１（図２Ａ）の各々について面積分を実行する必要が無いため、演算負荷が軽減される。
　上述の例では、体積要素５１を構成する複数のボクセル５０の各々に付与された磁化Ｍｊを用いて磁場を計算したが、粗視化後の体積要素５１の中心位置に、複数のボクセル５０に付与された磁気双極子を多重極展開して、多重極子に基づいて観測点Ｐｏの磁場を計算してもよい。
　双極子近似方式では、粗視化後の１つの体積要素５１に対応して１つの観測点Ｐｏが配置され、粗視化前の１つのボクセル５０が１つのソースとなる。このため、ｎ個のボクセル５０をまとめて１つの体積要素５１に粗視化すると、観測点Ｐｏの個数が１／ｎになり、ソースの個数は不変である。その結果、ソースと観測点との組み合わせの個数が１／ｎになることにより、磁場の演算時間がほぼｎに比例して減少する。
　［テーブル参照方式と双極子近似方式との比較］
　上述のテーブル参照方式と双極子近似方式とを比較すると、計算精度の点ではテーブル参照方式の方が優れており、演算負荷の点では双極子近似方式の方が優れている。なお、ソースとなるボクセル５０から遠く離れた領域の観測点の磁場を計算する際には、ボクセル５０内に分布する磁化を磁気双極子で近似しても十分な計算精度が得られる。
　従って、ソースとなるボクセル５０と観測点Ｐｏとの距離に基づいて、テーブル参照方式と双極子近似方式とのどちらを適用するか決定するとよい。例えば、ソースを中心としたある半径の球の内部の観測点Ｐｏの磁場の計算にはテーブル参照方式を適用し、この球の外側の観測点Ｐｏの磁場の計算には双極子近似方式を適用するとよい。
　［第３機能（逐次積分方式）］
　次に、図６を参照して、磁場演算部２３（図１）の第３機能（逐次積分方式）について説明する。
　図６は、ソースとなる体積要素５１Ｓと観測点Ｐｏとの位置関係を示す図である。観測点Ｐｏは、例えば体積要素５１の各々の中心に設定される。ソースとなる体積要素５１Ｓ内の位置ｒ’の磁化をＭ（ｒ’）で表す。観測点Ｐｏの位置をｒｉで表すと、ソースとなる体積要素５１Ｓに付与された磁化Ｍ（ｒ’）が観測点Ｐｏの位置に発生させる磁場Ｈ（ｒｉ）は、式（２）に基づいて計算することができる。ここで、式（２）の体積分は、ソースとなる体積要素５１Ｓ内の微小体積要素ｄＶ’について実行され、面積分は、ソースとなる体積要素５１Ｓの表面の微小領域ｄＳ’について実行される。
　体積要素５１Ｓ内の磁化Ｍ（ｒ’）が一定であると仮定できる場合、または体積要素５１Ｓ内における磁化Ｍ（ｒ’）の空間的な変化がほぼ無視できる場合、式（２）の右辺の体積分の項が０になるため、式（２）の右辺の面積分を実行することにより、磁場Ｈ（ｒｉ）を求めることができる。
　逐次積分方式においては、テーブル参照方式で利用した基準点Ｇｋの位置の磁場を平均化することなく、観測点Ｐｏの磁場を直接計算する。このため、テーブル参照方式に比べて、磁場の計算精度を高めることができる。特に、ソース点が観測点に近接しているときに、逐次積分方式を適用する効果が顕著である。
　また、テーブル参照方式では、ソースとなる体積要素５１Ｓ（図４）内の磁化Ｍｊが一定であるという条件の下で磁場が計算される。逐次積分方式は、体積要素５１内の磁化が空間的に変化している場合にも適用することができる。
　逐次積分方式では、粗視化後の１つの体積要素５１が１つのソースとなり、かつ１つの観測点となる。このため、ｎ個のボクセル５０をまとめて１つの体積要素５１に粗視化すると、ソース及び観測点の個数が共に１／ｎになる。その結果、ソースと観測点との組み合わせの個数が１／ｎ２になるため、磁場の演算時間がほぼｎ２に比例して減少する。
　次に、図７Ａを参照して、観測点Ｐｏに発生する磁場を計算する一例について説明する。
　図７Ａは、ソースＰｓと観測点Ｐｏ１、Ｐｏ２との位置関係を示す図である。ソースＰｓを中心とした半径ａ１の球の内側の第１領域Ｒ１に存在する観測点Ｐｏ１に発生する磁場の計算には、テーブル参照方式を適用する。これにより、計算精度を確保しつつ、演算負荷の軽減を図ることができる。第１領域Ｒ１の外側の第２領域Ｒ２に存在する観測点Ｐｏ２に発生する磁場の計算には、双極子近似方式を適用する。これにより、計算精度を確保しつつ、演算負荷をさらに軽減させることができる。
　次に、図７Ｂを参照して、観測点Ｐｏに発生する磁場を計算する他の例について説明する。
　図７Ｂは、ソースＰｓと、観測点Ｐｏ２、Ｐｏ３、Ｐｏ４との位置関係を示す図である。図７Ｂに示した例では、第１領域Ｒ１が、ソースＰｓを中心とする半径ａ２（ａ２＜ａ１）の球の内側の第３領域Ｒ３と外側の第４領域Ｒ４とに区分されている。第４領域Ｒ４内の観測点Ｐｏ４に発生する磁場の計算には、図７Ａの例と同様にテーブル参照方式を適用する。第３領域Ｒ３内の観測点Ｐｏ３に発生する磁場の計算には、逐次積分方式を適用する。これにより、計算精度をより高めることができる。
　実施例による磁場解析装置を用いることにより得られる効果を確認するために、実際に磁場解析を行った。次に、図８Ａ～図８Ｃを参照して、この解析結果について説明する。
　図８Ａ及び図８Ｂは、それぞれ解析対象の磁石電動機に用いられる永久磁石の平面図及び側面図である。ｚ軸を中心軸とする円柱の側面に沿って、円周方向に１０個の永久磁石６０が等間隔に配置されている。半径方向をｘ軸及びｙ軸と定義する。永久磁石６０の各々の平断面は円弧状であり、半径方向に磁化されている。円周方向に隣り合う２つの永久磁石６０の磁化方向は反対向きである。磁場解析において、永久磁石６０を、複数のボクセル５０（図２Ａ）に分割する。図８Ｂには、粗視化後の体積要素５１のｚ軸方向の境界を示しており、ｘ軸方向及びｙ軸方向の境界は示していない。
　永久磁石６０によって形成される磁場のｚ軸方向の変化は、ｘ軸方向及びｙ軸方向の変化に比べて緩やかである。このため、ｚ軸方向に連続する複数のボクセル５０（図２Ａ）をまとめて体積要素５１（図２Ｂ）とする粗視化を行う。また、ｚ方向の中央近傍において磁場の変化が相対的に緩やかであり、ｚ軸方向の両端近傍では磁場の変化が急峻になる。このため、１つの体積要素５１にまとめるボクセル５０の個数を、永久磁石６０のｚ軸方向の中心近傍において相対的に多くし、ｚ軸方向の両端近傍において相対的に少なくした。この磁場解析では、ボクセル５０の個数を２３４，９５２個とし、粗視化後の体積要素５１の個数を２６，０６４個とした。比較のために、粗視化しない比較例、すなわち、すべての体積要素５１の各々が１つのボクセル５０で構成される例についても、磁場解析を行った。
　また、磁場解析において、図７Ｂの第３領域Ｒ３と第４領域Ｒ４との境界となる球面の半径ａ２を、ボクセル５０の一辺の長さの４倍とし、第４領域Ｒ４と第２領域Ｒ２との境界となる球面の半径ａ１を、ボクセル５０の一辺の長さの１６倍とした。
　図８Ｃは、実施例による方法で磁場解析を行った結果、及び比較例による方法で磁場解析を行った結果を示すグラフである。横軸は、ｚ軸を中心とした回転方向の角度を単位「ラジアン」で表し、縦軸は磁場の強さを任意単位で表す。図８Ｃのグラフには、ｚ軸方向の中間における永久磁石６０の表面における磁場の強さを示している。グラフ中の実線は実施例による磁場解析結果を示し、丸記号は、比較例による磁場解析結果を示す。
　図８Ｃに示すように、実施例と比較例との磁場解析結果がよく一致している。これは、粗視化しても十分な計算精度が得られていることを意味する。比較例による方法で磁場解析を行ったときの計算時間は、実施例による方法で磁場解析を行ったときの計算時間の３３．９倍であった。粗視化することにより、演算負荷が軽減され、磁場解析時間の短縮化を図ることが可能であることが確認された。
　［双極子近似方式への多重極展開の適用］
　次に、図９～図１２を参照して、双極子近似方式を用いた磁場の計算に多重極展開を適用して計算時間の短縮を図る方法について説明する。
　図９は、多重極展開の方法を説明するために用いる座標系を示す図である。磁場の変化が相対的に緩やかな方向に一列に並ぶ複数の立方体のボクセル５０を粗視化して１つの体積要素５１が得られている。多重極展開する複数のボクセル５０の個数をＮで表す。粗視化された体積要素５１の中心の位置をｒｃで表す。位置ｒｃを多重極展開点Ｃということとする。多重極展開点Ｃを原点とし、磁場の変化が相対的に緩やかな方向をｚ軸方向とするｘｙｚ直交座標系を定義する。多重極展開点Ｃに多重極展開される対象となる磁気双極子を持つ複数のボクセル５０はｚ軸上に位置する。ｘｙｚ直交座標系において、ｊ番目のボクセル５０の位置はｒｊ－ｒｃで表され、観測点Ｐｏの位置はｒｉ－ｒｃで表される。ｊ番目のボクセル５０に付与された磁化をＭｊで表し、ｊ番目のボクセル５０の体積をΔＶｊで表す。ｚ軸と、多重極展開点Ｃから観測点Ｐｏに向かうベクトルｒｉ－ｒｃとのなす角度をθｉで表す。
　ｊ番目のボクセル５０と観測点Ｐｏとの間の距離ｒｉｊは、以下の式（８）で表される。
Figure JPOXMLDOC01-appb-M000008
 　ここで、ｒｉｃ、ｒｊｃは、それぞれ多重極展開点Ｃから観測点Ｐｏまで距離、及び多重極展開点Ｃからｊ番目のボクセル５０までの距離である。ＰＬ（ｃｏｓθｉ）は、ｃｏｓθｉを引数とするＬ次のルジャンドル多項式である。なお、式（８）は、J. D. JACKSON, CLASSICAL ELECTRODYNAMICS THIRD EDITION, p.103に説明されている。
　ベクトルｒｉ－ｒｃのｘ成分であるｘｉｃで式（８）を偏微分すると、以下の式（９）が得られる。
Figure JPOXMLDOC01-appb-M000009
 　ここで、ｒｉｃ＝ｒｉ－ｒｃ、ｒｊｃ＝ｒｊ－ｒｃである。ｘｊｃは、位置ベクトルｒｊｃのｘ成分である。
　ｊ番目のボクセル５０はｚ軸上に位置するため、ｘｊｃ＝ｙｊｃ＝０である。この条件を式（９）に適用すると、以下の式（１０）が導出される。
Figure JPOXMLDOC01-appb-M000010
 　ここで、ｇ１、ｇ２は以下の式で表される。
Figure JPOXMLDOC01-appb-M000011
 　また、ｚｊｃは、ベクトルｒｊｃのｚ成分である。
　式（１０）をｘｉｃで偏微分した後、－ｘｉｃで割ると、以下の式が得られる。
Figure JPOXMLDOC01-appb-M000012
 　ここで、ｆ１、ｆ２、ｆ３は、以下の式で表される。
Figure JPOXMLDOC01-appb-M000013
　磁気モーメントｍｊを持つｊ番目のボクセル５０が遠方の観測点Ｐｏに作る磁場Ｈ（ｒｉ）は、以下の式で近似することができる。
Figure JPOXMLDOC01-appb-M000014
 　ここで、Ｍｊ及びΔＶｊは、それぞれｊ番目のボクセル５０に付与された磁化及びｊ番目のボクセル５０の体積である。Ｎ個のボクセル５０が粗視化された体積要素５１が遠方の観測点Ｐｏに作る磁場は、式（１４）をＮ個のボクセル５０について足し合わせることにより求まる。式（１４）をＮ個のボクセル５０について足し合わせると、式（７）が得られる。
　式（１４）に式（１０）及び式（１２）を代入して整理すると、以下の式が得られる。
Figure JPOXMLDOC01-appb-M000015
 　ここで、中カッコ内の１段目、２段目、及び３段目の式は、それぞれｘ成分、ｙ成分、及びｚ成分を表している。ｍｊｚは、磁気モーメントｍｊのｚ成分である。
　粗視化前のＮ個のボクセル５０に対応する磁気モーメントｍｊが等しいと仮定する。このとき、ｍ１＝ｍ２＝・・・＝ｍＮ＝ｍと表すことができる。また、粗視化前のＮ個のボクセル５０の一辺の長さｄは等しいと仮定する。ｊ番目のボクセル５０の一辺の長さをｄｊで表すと、ｄ１＝ｄ２＝・・・＝ｄＮ＝ｄと表すことができる。観測点ＰｏへのＮ個のボクセルからの寄与を考えると、式（１５）から以下の式が導出される。
Figure JPOXMLDOC01-appb-M000016
 　ここで、ｍｚは、磁気モーメントｍのｚ成分である。係数ｓ１ Ｌ、ｓ２ Ｌ、ｓ３ Ｌは、以下の式で定義される。
Figure JPOXMLDOC01-appb-M000017
　粗視化の対象となる複数のボクセル５０がｚ軸方向にのみ配列しており、一辺の長さｄｊがすべてのボクセル５０において同一である場合、式（１７）に示す係数ｓ１ Ｌ、ｓ２ Ｌ、ｓ３ Ｌは、ボクセル５０を粗視化して体積要素５１を定義する段階で計算することが可能である。式（１６）には、位置ベクトルとしてｒｉｊが現れておらず、ｒｉｃのみが現れている。このため、粗視化後の体積要素５１が観測点Ｐｏに作る磁場は、多重極展開点Ｃの位置ｒｃと、観測点Ｐｏの位置ｒｉとを用いて計算することが可能になる。これは、観測点Ｐｏの磁場を計算するにあたり、粗視化前のＮ個のボクセル５０を個別に扱う必要がないことを意味する。従って、Ｎ個のボクセル５０の各々の位置ｒｊごとに磁場を計算する手法と比べて、計算時間を短縮することが可能である。
　従って、磁場の変化が相対的に緩やかな一方向（ｚ軸方向）に一列に配列した複数のボクセル５０の各々の磁気双極子を多重極展開点Ｃに多重極展開する手法を用いることにより、計算時間を短縮することが可能である。なお、計算時間の短縮を図るために、１つの体積要素５１に粗視化されるＮ個の複数のボクセル５０の一辺の長さｄｊが等しく、Ｎ個のボクセル５０にそれぞれ付与される磁化Ｍｊが等しいという条件を満たすようにシミュレーション条件を設定することが好ましい。
　係数ｓ１ Ｌ、ｓ２ Ｌ、ｓ３ Ｌの値は、実際にシミュレーション計算を行うときに用いるルジャンドル多項式の展開次数まで予め計算しておき、記憶装置２５（図１）に記憶させておくとよい。
　次に、多重極展開手法を用いることの効果を確認するために行った評価実験の結果について説明する。多重極展開の手法を用いた磁場計算と、個々のボクセル５０について双極子近似を適用した磁場計算とを実際に行って両者の精度及び計算時間を比較した。
　図１０Ａは、個々のボクセル５０について双極子近似を行い、観測点Ｐｏに発生する磁場Ｈ（ｒｉ）を計算により求める方法を説明する図である。１００個のボクセル５０がｚ軸方向に並んで配置されている。ボクセル５０の各々に磁化Ｍが付与されている。個々のボクセル５０が観測点Ｐｏに発生させる磁場を求め、１００個のボクセル５０について足し合わせることにより、磁場Ｈ（ｒｉ）を求めた。この計算によって求められた磁場をＨｏと表記することとする。
　図１０Ｂは、多重極展開手法を用い、観測点Ｐｏに発生する磁場Ｈ（ｒｉ）を計算により求める方法を説明する図である。１００個のボクセル５０（図１０Ａ）を粗視化して体積要素５１を定義し、体積要素５１の中心である多重極展開点Ｃに、１００個のボクセル５０の磁気モーメントを多重極展開した。立方体のボクセル５０の各々の一辺の長さをｄで表す。このとき、粗視化スケールは１００ｄになる。式（１７）を用いて係数ｓ１ Ｌ、ｓ２ Ｌ、ｓ３ Ｌを予め計算しておき、式（１６）を用いて観測点Ｐｏの磁場Ｈ（ｒｉ）を求めた。この計算によって求められた磁場をＨｃと表記することとする。磁場Ｈｃの計算は、ルジャンドル多項式の展開次数を０次、２次、４次とした３つの場合について行った。
　図１１は、磁場Ｈｏに対する磁場Ｈｃの誤差を、規格化距離の関数として示すグラフである。規格化距離は、多重極展開点Ｃから観測点Ｐｏまでの距離を粗視化スケール１００ｄで規格化したものである。誤差は、（Ｈｃ－Ｈｏ）／Ｈｏで定義される。図１１のグラフに示した細い実線、破線、及び太い実線は、それぞれルジャンドル多項式の展開次数を０次、２次、及び４次とした場合の誤差を示している。
　ルジャンドル多項式の展開次数を高くするほど、誤差が小さくなっていることがわかる。また、規格化距離が長くなるほど、誤差が小さくなっていることがわかる。図７Ａ及び図７Ｂに示したように、ソースＰｓから観測点Ｐｏまでの距離が短い範囲で逐次積分方式やテーブル参照方式を用いる場合は、ソースＰｓからある程度遠い範囲において双極子近似法式が適用される。このため、ルジャンドル多項式の展開次数を低くしても、十分高い精度が得られると考えられる。
　図１２は、多重極展開を行うことによる計算の高速化効果を示すグラフである。個々のボクセル５０について双極子近似を行う手法（図１０Ａ）を適用した場合の計算時間をＴ０で表し、多重極展開手法を適用した場合の計算時間をＴＣで表したとき、高速化効果をＴ０／ＴＣと定義する。例えば、計算時間ＴＣが計算時間Ｔ０の１／２になった場合、高速化効果は２倍である。
　ルジャンドル多項式の展開次数を低くするごとに、高速化効果が高くなっていることがわかる。ルジャンドル多項式の展開次数を２次までとした場合、約２５倍の高速化効果が得られている。また、ルジャンドル多項式の展開次数を４次までとした場合でも、約５倍の高速化効果が得られている。
　図１１及び図１２に示した評価実験により、多重極展開手法を適用することによって、所望の精度を維持しつつ、計算時間を短縮することが可能であることが確認された。ルジャンドル多項式の展開次数を低くすると計算時間は短くなるが、精度は低下する。求められる精度と、許容される計算時間に応じて、ルジャンドル多項式の展開次数を決めるとよい。
　図１１に示したように、規格化距離が長くなるとルジャンドル多項式の展開次数を低くしても高い精度が得られる。このため、規格化距離に応じてルジャンドル多項式の展開次数を変化させてもよい。規格化距離が長くなるに従って、ルジャンドル多項式の展開次数を低くするとよい。この手法により、求められる精度を維持しつつ、計算時間をさらに短くすることが可能になる。
　次に、上記実施例の変形例について説明する。
　上記実施例では、図２Ｂに示したｚ軸方向に関する磁場の変化が、ｘ軸方向及びｙ軸方向に関する磁場の変化より緩やかである場合を示した。この例は、磁化されている物体がｚ軸方向に長い棒状または筒状の形状を持つような場合に相当する。その他に、磁化されている物体がｘｙ面内に広がる板状の形状を持つような場合には、ｘ軸方向及びｙ軸方向に関する磁場の変化が、ｚ軸方向に関する磁場の変化より緩やかになる。このような場合には、体積要素５１のｘ軸方向及びｙ軸方向の寸法を、ｚ軸方向の寸法より大きくするとよい。
　上記実施例による解析手法は、ボクセル５０（図２Ａ）に内接する球形のビーズを配置し、ビーズに磁化（磁気モーメント）を付与し、ビーズの集合体からなる磁性材料の運動をシミュレーションする方法（磁気ビーズ法）にも適用することができる。磁気ビーズ法については、特開２０１５－１１１４０１号公報、特開２０１７－１９４８８４号公報に説明されている。
　式（７）に示した体積ΔＶｊを磁気ビーズの体積と考えると、ボクセル５０の体積は、磁気ビーズの充填率γｊを用いて、ΔＶｊ／γｊで表される。従って、本実施例による手法を磁気ビーズ法に適用する場合、式（７）のΔＶｊに代えて、ΔＶｊ／γｊを用いるとよい。
　図８Ａ～図８Ｃに示したシミュレーションでは、磁石電動機の永久磁石が磁化されている例を取り扱ったが、上記実施例による解析手法は、その他の電磁アクチュエータの磁場解析に利用することができる。また、対象物が永久磁石ではなく、鉄等の強磁性材料である場合にも、上記実施例による解析手法を適用することができる。この場合には、各ボクセル５０（図２Ａ）には、外部磁場によって生じる磁化を付与し、この磁化による磁場を求めればよい。
　上記実施例では、解析対象をボクセル５０（図２Ａ）に分割し、連続する複数のボクセル５０を粗視化して体積要素５１を形成したが、ボクセル５０を経ることなく、解析対象を直接体積要素５１に分割してもよい。
　上記実施例では、１つの体積要素５１（図２Ｂ）内に１つの観測点を設定したが、１つの体積要素５１内に複数の観測点を設定してもよいし、ボクセル５０ごとに観測点を設定してもよい。１つの体積要素５１内に複数の観測点を設定し、観測点に発生する磁場に基づいて磁化を生じさせると、１つの体積要素５１内に磁化の分布が生じる。テーブル参照方式を適用する場合には、１つの体積要素５１内の磁化の分布を平均化した磁化Ｍｊを用いて磁場を計算するとよい。双極子近似方式を適用する場合には、観測点の磁場から、ソースとなるボクセル５０（図５）ごとの磁化Ｍｊを求めるとよい。逐次積分方式を適用する場合には、ソースとなる体積要素５１内の均一ではない磁化Ｍ（ｒ’）に基づいて、観測点の磁場を計算すればよい。
　上述の実施例は例示であり、異なる例で示した構成の部分的な置換または組み合わせが可能であることは言うまでもない。複数の例の同様の構成による同様の作用効果については実施例ごとには逐次言及しない。さらに、本発明は上述の実施例に制限されるものではない。例えば、種々の変更、改良、組み合わせ等が可能なことは当業者に自明であろう。
２０　処理装置
２１　解析情報取得部
２２　磁化付与部
２３　磁場演算部
２４　出力制御部
２５　記憶装置
２６　データ保持部
２７　対応関係保持部
２８　入力装置
２９　出力装置
３０　磁場解析装置
５０　ボクセル
５１　体積要素
５１Ｓ　ソース点となる体積要素
６０　永久磁石

Claims (13)

  1. 　仮想空間を複数の体積要素に分割し、前記体積要素の各々に磁化を付与する磁化付与部と、
    　前記仮想空間の複数の観測点の各々について、前記観測点の周りの複数の前記体積要素に付与された磁化に基づいて、前記体積要素ごとに前記観測点に発生する磁場を演算し、複数の前記体積要素ごとの演算結果に基づいて、複数の前記観測点の各々に発生する磁場を求める磁場演算部と
    を有し、
    　前記磁化付与部は、少なくとも一部の前記体積要素について、第１方向に関する寸法を、他の第２方向に関する寸法より大きくする磁場解析装置。
  2. 　前記磁化付与部は、
    　前記仮想空間を、形状及び寸法が同一の複数のボクセルに分割し、
    　複数の前記ボクセルの各々に磁化を付与し、
    　前記第１方向に連続する複数の前記ボクセルをまとめて１つの前記体積要素とする粗視化を行い、
    　粗視化された前記体積要素に含まれる複数の前記ボクセルに付与されている磁化に基づいて、粗視化された前記体積要素に磁化を付与する請求項１に記載の磁場解析装置。
  3. 　前記体積要素に付与される磁化と、周りの複数の前記観測点に発生する磁場との対応関係を、前記体積要素と前記観測点との相対位置ごとに保持する対応関係保持部を、さらに有し、
    　前記磁場演算部は、前記対応関係を用いて、前記体積要素に付与されている磁化が前記観測点に発生させる磁場を求める第１機能を持つ請求項２に記載の磁場解析装置。
  4. 　前記対応関係保持部は、前記体積要素の大きさごとに前記対応関係を保持しており、
    　前記磁場演算部は、前記体積要素に付与されている磁化が前記観測点に発生させる磁場を計算する際に、磁場を発生させるソースとなる前記体積要素の大きさに対応する前記対応関係を用いて磁場を計算する請求項３に記載の磁場解析装置。
  5. 　前記磁場演算部は、複数の前記ボクセルに付与されている磁化をそれぞれ磁気双極子で代表して、前記磁気双極子に基づいて前記観測点の磁場を求める第２機能を有する請求項３乃至４のいずれか１項に記載の磁場解析装置。
  6. 　前記第２機能は、複数の前記ボクセルの前記磁気双極子によってそれぞれ前記観測点に発生する磁場を前記ボクセルごとに求め、前記ボクセルごとに求められた磁場を足し合わせる機能を含む請求項５に記載の磁場解析装置。
  7. 　前記第２機能は、一方向に一列に並ぶ複数の前記ボクセルの各々の前記磁気双極子を、１つの多重極展開点に多重極展開し、前記多重極展開点に対する前記観測点の相対位置に基づいて前記観測点の磁場を求める機能を含む請求項５に記載の磁場解析装置。
  8. 　１つの前記多重極展開点に多重極展開される複数の前記磁気双極子の大きさ及び向きは同一であり、１つの前記多重極展開点に多重極展開される複数の前記磁気双極子に対応する前記ボクセルの寸法は同一である請求項７に記載の磁場解析装置。
  9. 　前記磁場演算部は、磁場を発生させるソースとなる前記体積要素から前記観測点までの距離が第１距離以下の少なくとも一部の前記観測点に発生する磁場を、前記第１機能を用いて求め、前記第１距離より長い前記観測点に発生する磁場を、前記第２機能を用いて求める請求項５乃至８のいずれか１項に記載の磁場解析装置。
  10. 　前記磁場演算部は、磁場を発生させるソースとなる前記体積要素に付与された磁化を含む物理量を、前記体積要素の表面について面積分することにより、前記体積要素の周りの前記観測点に発生する磁場を求める第３機能を有する請求項９に記載の磁場解析装置。
  11. 　前記磁場演算部は、磁場を発生させるソースとなる前記体積要素から前記観測点までの距離が前記第１距離より短い第２距離以下の前記観測点に発生する磁場を、前記第３機能を用いて求める請求項１０に記載の磁場解析装置。
  12. 　仮想空間を複数の体積要素に分割し、前記体積要素の各々に磁化を付与し、
    　前記体積要素に付与された磁化が前記仮想空間内の複数の観測点に発生させる磁場を前記体積要素ごとに計算し、複数の前記体積要素ごとの計算結果に基づいて、複数の前記観測点の各々に発生する磁場を求め、
    　前記仮想空間を前記体積要素に分割するにあたり、少なくとも一部の前記体積要素について、前記仮想空間内の第１方向に関する寸法を、他の第２方向に関する寸法より大きくする磁場解析方法。
  13. 　仮想空間を複数の体積要素に分割し、前記体積要素の各々に磁化を付与する機能と、
    　前記体積要素に付与された磁化が前記仮想空間内の複数の観測点に発生させる磁場を前記体積要素ごとに計算し、複数の前記体積要素ごとの計算結果に基づいて、複数の前記観測点の各々に発生する磁場を求める機能と、
    　前記仮想空間を前記体積要素に分割するにあたり、少なくとも一部の前記体積要素について、前記仮想空間内の第１方向に関する寸法を、他の第２方向に関する寸法より大きくする機能と
    をコンピュータに実現させるためのプログラム。
PCT/JP2018/041831 2018-03-07 2018-11-12 磁場解析装置、解析方法、及びプログラム WO2019171660A1 (ja)

Priority Applications (2)

Application Number Priority Date Filing Date Title
JP2020504776A JP7121111B2 (ja) 2018-03-07 2018-11-12 磁場解析装置、解析方法、及びプログラム
US17/012,682 US11703553B2 (en) 2018-03-07 2020-09-04 Magnetic field analysis device, analysis method, and computer readable medium storing program

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2018040400 2018-03-07
JP2018-040400 2018-03-07

Related Child Applications (1)

Application Number Title Priority Date Filing Date
US17/012,682 Continuation US11703553B2 (en) 2018-03-07 2020-09-04 Magnetic field analysis device, analysis method, and computer readable medium storing program

Publications (1)

Publication Number Publication Date
WO2019171660A1 true WO2019171660A1 (ja) 2019-09-12

Family

ID=67846017

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2018/041831 WO2019171660A1 (ja) 2018-03-07 2018-11-12 磁場解析装置、解析方法、及びプログラム

Country Status (3)

Country Link
US (1) US11703553B2 (ja)
JP (1) JP7121111B2 (ja)
WO (1) WO2019171660A1 (ja)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2021131766A (ja) * 2020-02-20 2021-09-09 住友重機械工業株式会社 磁場解析方法、磁場解析装置、及びプログラム

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2021110997A (ja) * 2020-01-07 2021-08-02 住友重機械工業株式会社 シミュレーション方法、シミュレーション装置、及びプログラム

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2003006181A (ja) * 2001-06-26 2003-01-10 Fujitsu Ltd 電磁波解析装置および電磁波解析プログラム

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP5052985B2 (ja) 2007-07-31 2012-10-17 住友重機械工業株式会社 分子シミュレーション方法、分子シミュレーション装置、分子シミュレーションプログラム、及び該プログラムを記録した記録媒体
JP6249912B2 (ja) 2013-11-01 2017-12-20 住友重機械工業株式会社 解析装置
US10222441B2 (en) * 2016-04-03 2019-03-05 Q Bio, Inc. Tensor field mapping
JP6618417B2 (ja) 2016-04-22 2019-12-11 住友重機械工業株式会社 解析装置および解析方法

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2003006181A (ja) * 2001-06-26 2003-01-10 Fujitsu Ltd 電磁波解析装置および電磁波解析プログラム

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2021131766A (ja) * 2020-02-20 2021-09-09 住友重機械工業株式会社 磁場解析方法、磁場解析装置、及びプログラム
JP7264841B2 (ja) 2020-02-20 2023-04-25 住友重機械工業株式会社 磁場解析方法、磁場解析装置、及びプログラム

Also Published As

Publication number Publication date
US20200400757A1 (en) 2020-12-24
US11703553B2 (en) 2023-07-18
JP7121111B2 (ja) 2022-08-17
JPWO2019171660A1 (ja) 2021-03-11

Similar Documents

Publication Publication Date Title
JP6249912B2 (ja) 解析装置
WO2019171660A1 (ja) 磁場解析装置、解析方法、及びプログラム
Michoski et al. Foundations of the blended isogeometric discontinuous Galerkin (BIDG) method
JP6618417B2 (ja) 解析装置および解析方法
JP5839473B2 (ja) 解析装置
JP2005100067A (ja) マイクロ磁化解析プログラムおよび解析装置
JP5546439B2 (ja) 解析装置および解析方法
JP5570307B2 (ja) 解析装置および解析方法
Gusenbauer et al. Automated meshing of electron backscatter diffraction data and application to finite element micromagnetics
JP5546440B2 (ja) 解析装置および解析方法
Weeber et al. Accelerating the calculation of dipolar interactions in particle based simulations with open boundary conditions by means of the P2NFFT method
Simona Numerical methods for the simulation of particle motion in electromagnetic fields
JP5236599B2 (ja) 磁場解析装置および磁場解析装置の動作方法
JP5236601B2 (ja) 磁場解析装置および磁場解析装置の動作方法
JP2011244623A (ja) 連成解析システム、解析システム、および連成解析方法
JP5460181B2 (ja) 磁場解析装置および磁場解析装置の動作方法
JP2022112214A (ja) 解析方法、解析装置、及びプログラム
CN112069742A (zh) 使显式数值方案稳定化
JP5478154B2 (ja) 磁場解析装置および磁場解析装置の動作方法
Shivaram et al. A finite element mesh generation for complex geometry
Cheviakov Construction of exact plasma equilibrium solutions with different geometries
Farrahi et al. Efficiency of hysteresis rods in small spacecraft attitude stabilization
Sommer et al. Interactive high-resolution simulation of granular material
Rouse Three dimensional computer modeling of electron optical systems
JP7264841B2 (ja) 磁場解析方法、磁場解析装置、及びプログラム

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 18908511

Country of ref document: EP

Kind code of ref document: A1

WWE Wipo information: entry into national phase

Ref document number: 2020504776

Country of ref document: JP

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 18908511

Country of ref document: EP

Kind code of ref document: A1