WO2017183172A1

WO2017183172A1 - Computer and calculation method

Info

Publication number: WO2017183172A1
Application number: PCT/JP2016/062715
Authority: WO
Inventors: 辰也戸丸
Original assignee: 株式会社日立製作所
Priority date: 2016-04-22
Filing date: 2016-04-22
Publication date: 2017-10-26
Also published as: JPWO2017183172A1; JP6524340B2

Abstract

According to the present invention, the accuracy of solutions obtained by a local field response method, in which spins serving as variables are caused to be oriented in the direction of a local effective magnetic field, is improved by using a discrete time axis and setting a spin response function with respect to an effective magnetic field so as to depend on two adjacent times. Specifically, N variables sj^z have the range -1 ≤ sj^z ≤ 1, and a problem is set up using a coefficient gj, which represents a local item, and another coefficient Jkj, which represents variable-variable interaction. Then calculation is performed for time t which is discrete, namely t = t0 (t0 = 0) to tm (tm = τ) obtained by dividing a time period by m. For each time ti, variables Bj^z0 (ti), Bj^z(ti), and sj^z(ti) are sequentially determined by performing the following calculations: Bj^z0(ti) = Σk(≠j)Jkjsk^z(t_i-1) + gj Bj^z(ti) = (1-u)Bj^z0(ti) + uBj^z(t_i-1) Beff,j^z(ti) = Bj^z(ti)∙ti/τ sj^z(ti) = f(Beff, j^z(ti), ti) (where: Bj^z(t0) = 0; sj^z(t0) = 0; u is a parameter satisfying the equation 0 ≤ u ≤ 1; and f is a function defined such that sj^z(ti) has the range -1 ≤ sj^z(ti) ≤ 1). For each successive time step, that is, as the time t is incrementally increased from t0 to tm, the variable sj^Z is caused to approach -1 or 1. At the end of this calculation, a solution is found by assuming that sj^zfd = -1 if sj^z < 0, or assuming that sj^zfd = 1 if sj^z > 0.

Description

Calculator and calculation method

The present invention relates to a computer that enables high-speed computation for inverse problems and combinatorial optimization problems that require exhaustive search, and a calculation method using the computer.

現代 Today, huge data processing is required as represented by the words “big data”. Data processing methods are expected to develop in the future, including social scientific and economic meanings, but the processing itself is the responsibility of the computer. The computer gives initial values and calculates based on the algorithm, and is good at forward calculation, but it is an inverse problem that estimates the initial value from the result and a combinatorial optimization problem that selects the optimal solution from many possibilities Is generally not good at all because it requires an exhaustive search in the worst case. However, for big data that requires various processing, processing such as exhaustive search is indispensable.

In recent years, a technique called quantum annealing, also known as adiabatic quantum computation, has been attracting attention (Non-Patent Documents 1 and 2). In this method, a problem is set so that the ground state of a certain physical system becomes a solution, and the solution is obtained through finding the ground state. Let H ^ _{p be} the Hamiltonian of the physical system that sets the problem. However, the Hamiltonian at the start of the calculation is not H ^ _p , but another Hamiltonian H ^ ₀ that is easy to prepare for the ground state. Next, take enough time to move the Hamiltonian from H ^ ₀ to H ^ _p . If enough time is taken, the system stays in the ground state, and finally the ground state (solution state) of the Hamiltonian H ^ _p is obtained. This is the principle of quantum annealing.

The ground state search method of physical system called Ising spin glass can cope with the problem called NP difficulty (Non-patent Document 3). Among combinatorial optimization problems, problems with high difficulty belong to NP difficulty. Furthermore, problems classified as P in computational complexity theory and problems classified as NP can all be reduced to NP-hard problems. Therefore, if quantum annealing is applied in the Ising spin glass system, almost all combinatorial optimization problems can be solved, greatly contributing to the processing of big data.

Another reason why quantum annealing is attracting attention is its robustness against decoherence. In quantum computers, quantum coherence had to be maintained over the computation time. On the other hand, in quantum annealing, a correct answer can be obtained if the ground state is maintained. It is not always necessary to maintain quantum coherence. Considering that it is difficult to construct a pure quantum system at the current technical level, and thus it is difficult to maintain quantum coherence over the calculation time, the reason why quantum annealing is attracting attention can be understood. However, quantum annealing also has drawbacks. The reason why quantum annealing can be realized is limited to the superconducting magnetic flux qubit system at present (Patent Document 1, Non-Patent Document 4) and requires a cryogenic cooling device. The necessity of cryogenic temperature is a challenge for realizing a practical computer.

A method devised to solve this problem is the local field response method described below (Patent Document 2, Non-Patent Document 5). First, reconsidering quantum annealing. The concept of annealing (annealing) originally exists regardless of quantum or classic, and quantum annealing uses quantum properties to improve the performance of classical annealing. This is why the quantum coherence does not necessarily have to be maintained over the calculation time in the quantum annealing, and the ground state should be maintained. There can be a methodology different from quantum annealing if the concept of annealing is applied regardless of quantum or classic. The above-mentioned local field response method was invented from that viewpoint. In this method, as in quantum annealing, a transverse magnetic field is applied to a spin system as a computing unit at time t = t ₀ , the magnetic field is gradually reduced, and a solution is obtained at time t = τ. In this method, the computing unit itself is classic, and quantum mechanical information is added to the spin response function to the magnetic field. Since this method is operated by a classic machine, it can be operated at room temperature and solves the problem of quantum annealing that requires extremely low temperature, but it is necessary to determine in advance a response function that reflects the quantum effect. In Patent Document 2 and Non-Patent Document 5, a response function including an average quantum effect is determined from the result of solving an empirical or similar problem. Although it contains an average quantum effect, the performance is improved compared to the purely classical method. However, in order to further improve the performance, it is necessary to make the response function dependent on individual problems and states.

Special table 2009-524857 International Publication WO2015 / 118639

As described above, quantum annealing requires a cryogenic cooling device in order to use a superconducting flux qubit. In addition, although the local field response method operates at room temperature, the performance improvement is limited because the quantum effect that brings about the performance improvement is average. SUMMARY OF THE INVENTION An object of the present invention is to provide a computer capable of operating at room temperature with sufficient performance for difficult problems that require exhaustive search.

One aspect of the present invention is that in a local field response method in which a spin as a variable is made to respond to a local effective magnetic field, the time axis is made discrete and the response function of the spin to the effective magnetic field is made to depend on two consecutive times. Is to make the time evolution similar to the quantum mechanical time evolution. More specifically, it is as follows.

N spin variables s _j ^z (j = 1, 2,…, N) take a range of -1 ≤ s _j ^z ≤ 1, and a coefficient g _j representing a local term and a coefficient J _kj representing an interaction between variables (k, j = 1, 2,…, N)
The time is divided into m and discretely calculated from t = t ₀ (t ₀ = 0) to t _m (t _m = τ)
At each time t _i (i = 1, 2, .., N), the variables B _j ^z0 (t _i ) and B _j ^z (t _i ) related to the effective magnetic field and the spin variable s _j ^z (t _i ) are The initial values at time t ₀ are B _j ^z (t ₀ ) = 0 and s _j ^z (t ₀ ) = 0, and B _j ^z0 (t _i ) is B _j ^z0 (t _i ). = (Σ _{k (≠ j)} J _kj s _k ^z (t _i−1 ) + g _j ), and B _j ^z (t _i ) is expressed as B _j ^z ( t _i ) = (1−u) B _j ^z0 (t _i ) + uB _j ^z (t _i−1 ), and B _j ^z (t _i ) is multiplied by factor t _i / τ to give B _{eff, j} ^z ( t _i ) = B _j ^z (t _i ) ・ t _i / τ
The s _j ^z (t _i ) is defined by s _j ^z (t _i ) = f (B _{eff, j} ^z (t _i ), t _i ) using the function f, and the function f is s _j ^z ( t _i ) is defined such that the range of −1 ≦ s _j ^z (t _i ) ≦ 1,
Close time step to -1 or 1 the variable s _j ^z as proceeding from t = t ₀ to t = t _m, finally s _j ^z <0 if _{^{_{s j zfd = -1, s j}}} z> 0 Then, it is characterized by s _j ^zfd = 1.

Another aspect of the present invention is one aspect of a computer including an input device, an output device, a storage device, a general arithmetic device, and a local field response arithmetic device. Here, N variables s _j ^z (j = 1, 2, ..., N) take a range of -1 ≤ s _j ^z ≤ 1, and a coefficient g _j representing a local term and a coefficient representing an interaction between variables The task is set by J _kj (k, j = 1, 2,…, N). The local field response calculation device divides the time into m and calculates discretely from t = t ₀ (t ₀ = 0) to t _m (t _m = τ), and each time t _i (i = 1, 2, .., N), the variables B _j ^z0 (t _i ), B _j ^z (t _i ), and s _j ^z (t _i ) are determined in order, and the initial value at time t ₀ is B _j ^z ( t ₀ ) = 0 and s _j ^z (t ₀ ) = 0, and B _j ^z0 (t _i ) is B _j ^z0 (t _i ) = (Σ _{k (≠ j)} J _kj s _k ^z (t _i-1 ) + g _j ) and B _j ^z (t _i ) is B _j ^z (t _i ) = (1-u) Bj ^z0 (t _i ) + uB _j ^z (t _i-1 ), B _j ^z (t _i ) is multiplied by factor t _i / τ to give B _{eff, j} ^z (t _i ) = B _j ^z (t _i ) · t _i / τ, and s _j ^z (t _i ) is defined by s _j ^z (t _i ) = f (B _{eff, j} ^z (t _i ), t _i ) using the function f, and the function f is s _j ^z (t _i range of) is defined to be _{^{_{-1 ≦ s j z (t i}}} ) ≦ 1. Then, as the time step is advanced from t = t ₀ to t = t _m , the variable s _j ^z approaches −1 or 1, and in a general arithmetic unit, if s _j ^z (t _i ) <0 at each time t _i If s _j ^zd (t _i ) = −1, s _j ^z (t _i )> 0 then s _j ^zd (t _i ) = 1 and if s _j ^z (t _i ) = 0, then s _j ^zd (t _i ) = 0, H _p (t _k ) =-Σ _{k> j} J _kj s _k ^zd (t _i ) s _j ^zd (t _i ) ^-Σ _j g _j s _j ^zd (t _i ) is calculated at each time t _i Then, ^let s _j ^zfd = s _j ^zd (t _{i ′} ) at time t _{i ′} when H _p (t _i ) becomes the minimum value be the final solution.

Another aspect of the present invention is a calculation method using a computer that includes a calculation unit, a storage unit, and a control unit, and performs calculations while exchanging data between the storage unit and the calculation unit under the control of the control unit. As one aspect. In this method, N variables s _j ^z (j = 1, 2, ..., N) take a range of -1 ≤ s _j ^z ≤ 1, and represent a coefficient g _j representing a local term and an interaction between variables. The task is set by the coefficient J _kj (k, j = 1, 2,…, N). In the calculation unit, the time is divided into m and discretely calculated from t = t ₀ (t ₀ = 0) to t _m (t _m = τ), and each time t _i (i = 1, 2,. , N), the variables B _j ^z0 (t _i ), B _j ^z (t _i ), and s _j ^z (t _i ) are determined in order, and the initial value at time t ₀ is B _j ^z (t ₀ ) = 0 and s _j ^z (t ₀ ) = 0, and B _j ^z0 (t _i ) is B _j ^z0 (t _i ) = (Σ _{k (≠ j)} J _kj s _k ^z (t _i− 1) + g _j ), and B _j ^z (t _i ) is B _j ^z (t _i ) = (1-u) B _j ^z0 (t _i ) + uB _j ^z ( t _{i- 1)} and then, the _{^{_{B j z (t i) B}}} eff is multiplied by the factor t _i / _τ, and _{^{_{j z (t i) = B}}} j z (t i) · t i / τ, s j z (t _i ) is defined by s _j ^z (t _i ) = f (B _{eff, j} ^z (t _i ), t _i ) using the function f, and the function f of s _j ^z (t _i ) The range is defined so that −1 ≦ s _j ^z (t _i ) ≦ 1. The closer the time step to -1 or 1 the variable s _j ^z as proceeding from t = t ₀ to t = t _m, finally s _j ^z <0 if _{^{_{s j zfd = -1, s j}}} z> If 0, set s _j ^zfd = 1 to determine the solution.

Another aspect of the present invention is one aspect of a calculation program itself, which is software stored in a storage unit, or a storage medium storing the same in order to cause the calculation unit to execute the above calculation method.

According to the present invention, a practical computer that can solve difficult problems with high accuracy is realized.

It is the conceptual diagram which showed the principle of this invention typically. 5 is a flowchart illustrating an example of an algorithm according to the first embodiment. Is a graph showing specific values of the response function r _b. It is the flowchart which showed the example of the algorithm which concerns on Example 2 with the flowchart. It is the flowchart which showed the example of the algorithm which concerns on Example 3 with the flowchart. It is the flowchart which showed the example of the algorithm which concerns on Example 4 with the flowchart. 10 is a flowchart illustrating an example of an algorithm according to a fifth embodiment. 10 is a flowchart illustrating an example of an algorithm according to a sixth embodiment. FIG. 10 is a block diagram illustrating an example of a computer configuration according to a seventh embodiment. It is the block diagram which showed an example regarding the detail of the local field response calculating apparatus part in the computer which concerns on Example 7. FIG.

Hereinafter, various embodiments of the present invention will be described with reference to the drawings, including the principle of calculation. However, the present invention is not construed as being limited to the description of the embodiments below. Those skilled in the art will readily understand that the specific configuration can be changed without departing from the spirit or the spirit of the present invention.

In the structure of the invention described below, the same portions or portions having similar functions are denoted by the same reference numerals in different drawings, and redundant description may be omitted.

In this specification and the like, notations such as “first”, “second”, and “third” are attached to identify the constituent elements, and do not necessarily limit the number or order. In addition, a number for identifying a component is used for each context, and a number used in one context does not necessarily indicate the same configuration in another context. Further, it does not preclude that a component identified by a certain number also functions as a component identified by another number.

The position, size, shape, range, etc. of each component shown in the drawings and the like may not represent the actual position, size, shape, range, etc. in order to facilitate understanding of the invention. For this reason, the present invention is not necessarily limited to the position, size, shape, range, and the like disclosed in the drawings and the like.

Example 1 starts with a quantum mechanical description and describes the principle that forms the basis of this example through the transition to the classical form.

FIG. 1 schematically shows the principle of this embodiment. The basic framework is the same as the local field response method described in Patent Document 2 and Non-Patent Document 5. Apply a transverse magnetic field at t = 0 to align the spins in one direction. After that, the transverse magnetic field is slowly reduced to the Hamiltonian of the problem setting at t = τ. The spin evolves in time in response to the local effective magnetic field applied at each time.

問題 The Hamiltonian of the problem setting and the Hamiltonian at t = 0

And the Hamiltonian at time t is

And τ is the calculation time. From the analogy of a one-spin system, the effective magnetic field received by the spin at site j is given by B ^ _{eff, j} = -∂H ^ / ∂σ ^ _j .

The local field response method of this embodiment is to operate on a classical machine by regarding <σ ^ _j > taking an expected value as a spin variable. As is clear from equations (1) and (2), <σ ^ _j > and consist only of x and z components. So the response function r _b a (t)

The spin direction is determined based on this response function. Spin system Motomari only effective magnetic field of each spin responses each site if classic, response function becomes r _b (t) = 1. However, quantum mechanics has non-local correlation and generally r _b (t) ≠ 1. As already mentioned, equation (5) has shifted to the classical equation by taking the expected value, but the quantum effect is taken in through r _b (t) ≠ 1. The value of r _b (t) is obtained by empirical or quantum mechanical precalculation of similar problems. The quantum effect here is average to be based on empirical or similar problems. The local field response method may be operated with r _b (t) = 1 without incorporating the quantum effect. The local field response method itself works without including the quantum effect.

Four variables <σ ^ _j ^z (t _i )>, <σ ^ _j ^x (t _i )>, , takes the expected value and is a classical quantity. Therefore, the notation of quantum mechanics is changed to the notation of classical physics. That is, <σ ^ _j ^x (t _i )> → s _j ^x (t _i ) <σ ^ _j ^z (t _i )> → s _j ^z (t _i ), → B _{eff, j} ^x (t _i ) → B _{eff, j} ^z (t _i ) With this change in notation, equation (5) becomes

It becomes.

Time evolution of the local field response method discretely performed, B _eff at time _{_{^{_{t i, j z (t i}}}} ) is determined from _{^{_{s k z (t i-1}}} ) at time t _i-1 according to equation (4). S _j ^z at time t _i (t _i) decide B _eff at time t _i in accordance with Equation _{(6), j} ^z a (t _i) to the original. Repeat this procedure. FIG. 2 summarizes this as a flowchart.

101 in FIG. 2 represents the setting of the initial value as the starting point of the algorithm. In 102a, B _{eff, j} ^z (t _i ) is obtained using s _k ^z (t _i−1 ) at time t _i−1 . The transverse magnetic field strength B _{eff, j} ^x (t _i ) is determined depending on the time. 103, B _{eff, j} ^z (t _i ) / B _{eff, j} ^x (t _i ) and response function r _b (t), corresponding to the spin direction, tanθ = r _b (t) ・ B _{eff, j} ^z ( t _i ) / B _{eff, j} ^x (t _i ) is obtained, and s _j ^z (t _i ) = r _s (t) · sin θ is obtained using the parameter r _s (t) representing the magnitude of the spin. r _s (t) is a quantity defined by r _s (t) ² = s _j ^x (t _i ) ² + s _j ^z (t _i ) ² and is empirically or precomputed like r _b (t) Decide on. Specific examples will be described in the next two paragraphs. 102a and 103 in FIG. 2 are one set of repeated calculation. Repeat this set until _{t = τ (= t m)} , s in _{^{_{t = τ j z (t m}}} )> 0 if _{^{_{s j zfd = 1, s j}}} z (t m) <0 if s _j ^zfd = The solution s _j ^zfd is obtained as −1 (processing 201 and 202).

Step 103 _{^{_{s j z (t i) =}}} f (B eff, j z (t i), t i) using the function f can be generally Writing _{^{as, f (B eff, j z}} ( t _i ), t _i ) = r _s (t) · sin {arctan (r _b (t) · B _{eff, j} ^z (t _k ) / B _{eff, j} ^x (t _k ))}. Since r _s (t) is a parameter representing the magnitude of the spin, 0 ≦ r _s (t) ≦ 1. Further, −1 ≦ s _j ^z (t _i ) ≦ 1. Here, r _b (t) = 1 and r _s (t) = 1 are purely classical.

FIG. 3 shows an example of pre-calculation of the response function r _b (t). This is a case where J _ij and g _j are determined by a uniform random number of [−5, 5] in an 8-bit system. It is the result of 100 problems and consists of 800 points (100 problems x 8 bits). Here, B _zx ≡ / and s _zx ≡ <σ ^ _j ^z > / <σ ^ _j ^x >. The point is a value determined strictly quantum mechanically. Reflecting the non-local correlation of quantum mechanics, the response function varies greatly. The circle is an average obtained by dividing the horizontal axis into 40 parts. Response functions are averaged is smooth B _zx dependent. If smooth, it can be described with several parameters. Non-Patent Document 5 describes a method for describing a smooth response function using four parameters, and the solid line r _b ⁰ (t) in FIG. 3 is obtained by this method. Another parameter r _s (t) is also obtained from the same four parameters.

FIG. 3 shows an example of the response function, and non-patent document 5 is cited to refer to the determination method of r _b ⁰ (t) and r _s (t). However, the response function r _b ⁰ (t) and r _s determination method (t) can also be empirically determined with there can be a variety of ways without being limited thereto. In the following embodiments, r _b ⁰ (t) is not limited to the solid line in FIG. 3 and represents an average response function.

Note that, similarly to r _b (t), r _s (t) can also be operated with r _s (t) = 1. Since the range of r _b (t) is −∞ <r _b (t) <∞, 0 ≦ r _s (t) ≦ 1, so the effect on the final solution by setting r _s (t) = 1 Is smaller than the effect of r _b (t) = 1. Therefore, setting r _s (t) = 1 is an effective means when prior information for determining r _s (t) is insufficient.

In the algorithm of FIG. 2, the calculation is repeated until t = τ (= t _m ), and the final solution s _j ^zfd is determined based on the sign of s _j ^z (t _m ). However, an optimal solution is not always obtained at t = τ. After an optimal solution is obtained at t <τ, the accuracy of the solution may decrease.

FIG. 4 shows a corresponding algorithm. As shown in the procedure 300 of FIG. 4, the energy is calculated at each time t = t _i , and the spin array to which the lowest energy is given through all times may be used as the final solution.

In Step 300 by the code decision was carried out s _j ^z at each time t = t _i in _{^{_{s j z (t i) (}}} t i)> 0 if _{^{_{s j zd (t i) =}}} 1, s j z (t i ) <0, s _j ^zd (t _i ) = − 1. When s _j ^z (t _i ) = 0, s _j ^zd (t _i ) = 0. (Processing 301 and 302). Next, an energy value H _p (t _i ) with respect to the value of s _j ^zd (t _i ) is obtained (step 303). Further, H _p (t _i ) is compared with H _p (t _i−1 ). _{_{H p (t i) <H}} p (t i-1) if _{_{^{t ≦ t i s j zd (}}} t i) the s _j range in H _p (t _i) since the lowest energy t = t _i in ^Save as ^zfd . If this process is ^repeated until t = t _m , an optimal solution remains in s _j ^zfd . Incidentally, the H _p (t ₀₎ = ₀ as an initial value.

Fig. 3 shows an example in which the response function is calculated quantum mechanically. Reflecting the non-local correlation of quantum mechanics, the response function varies greatly. In order to improve the performance of the local field response method, it is necessary to incorporate this variation into the response function. In this embodiment, the method will be described.

An important time zone in the calculation of the local field response method is the vicinity where the sign of the spin changes. In its vicinity is _{^{s j z (t) ≒ 0}} , B eff in conjunction _therewith, the _j ^z (t) ≒ 0. “≈” means substantially equal. The response function in Fig. 3 has been rewritten as equation (5).

It is. is in the denominator, so if ≈0, r _b (t) fluctuates greatly. Since the local field response method evolves according to the response function, the accuracy of the solution depends on whether r _b (t) fluctuates greatly or not.

The time at which the sign of the spin changes is in a linearly coupled state in terms of quantum mechanics, and the +1 state and the -1 state overlap at 50:50. In terms of band theory, this time is the band rebound point. Therefore, this time is also important in terms of quantum mechanics, and it is necessary to correctly incorporate the behavior in the vicinity of s _j ^z (t) ≈0 and B _{eff, j} ^z (t) ≈0 into the algorithm.

The linearly coupled state is one of the features unique to quantum mechanics, so it cannot be introduced into a classic machine as it is. Therefore, as described below, the effective magnetic field is determined from the spin value at two times, thereby causing a behavior similar to linear combination.

The effective magnetic field B _eff in Example 1 at time t _{_i,} was determined by _j ^z (t _i) the time t spin values in _{_{^{i-1 s j z (t}}} i-1). That is,

Is the fundamental factor in the effective magnetic field B _{eff, j} ^z (t _i ). Factors be taken into account at time t _i-1 of the previous one as well the time t _i in this embodiment

The effective magnetic field is determined based on Here, u is appropriately determined so that accuracy is high when 0 ≦ u ≦ 1. A typical value is u≈0.1. If we describe the effective magnetic field including the transverse magnetic field and its schedule,

It is. Around the time when the sign of the spin is reversed (the band repulsion point in quantum mechanics), the sign of the effective magnetic field is also reversed. Therefore, in a new effective magnetic field according to Equation (9), B _j ^z0 (t _i−1 ) and B _j ^z0 (t _i ) cancel each other, and B _j ^z (t _i ) ^≈0 . The spin value is calculated using this new effective magnetic field and the averaged response function r _b ^0.

Determined by Thus, if B _j ^z (t _i ) ≈0, s _j ^z (t) ≈0.

In quantum mechanics, the spin state is a linearly coupled state. On the other hand, in the handling here, the effective magnetic field defined by Equation (9) is a linear combination of the effective magnetic fields at two times. Therefore, they do not do the same thing mathematically. However, handling of equation (9) makes it possible to reproduce behavior similar to quantum mechanical spin on a classical machine. This improves the accuracy of the local field response solution.

In this treatment in which the effective magnetic field depends on two times, an average response function r _b ⁰ is used as a response function as shown in Equation (11). On the other hand, the response function corresponding to Equation (7) is obtained by using the effective magnetic field B _{eff, j} ^z0 (t _i ) = (t _i / τ) B _j ^z0 (t _i ) determined only by time t _i.

It becomes. This r _b varies due to the difference in the definitions of Equation (11) and Equation (12). The response function based on Equation (12) qualitatively reproduces the dispersed response function illustrated in FIG.

FIG. 5 shows a flowchart based on the above principle. The difference from FIG. 4 is that the procedure 102a is changed to the procedure 102b. This change is based on _Eqs. (8), (9), and (10) _.The effective magnetic field B _{eff, j} ^z (t _i ) is improved, and the variation of the response function r _b (t _i ) is reproduced. The Further, in FIG. 4, the response function in the procedure 103 is described by general r _b (t), but in FIG. 5, r _b ⁰ (t) is set according to the equation (11).

In the present embodiment, in order to reproduce a phenomenon in the vicinity of s _j ^z (t) ≈0, B _{eff, j} ^z (t) _≈0, the handling of Expression (9) is introduced. This treatment has quantum mechanical significance outside of this area. As a simple example, consider a two-spin system with J ₁₂ > 0 and g ₁ = g ₂ = 0. If _{^{_{s 1 z (t i)>}}} 0, s 2 z (t i) < assumed to be 0. In this case equation (8) by and be determined to the effective magnetic field the effective magnetic field only in their opponent spins is determined both spin sign inversion _{^{_{s 1 z (t i + 1}}} ) <0, s 2 z (t _{i + 1} )> 0. After that, both spins repeat inversion.

However, in an actual quantum mechanical spin system, the subsequent spin state is determined depending not only on the spin state of the partner but also on its own spin state. Equation (9) has this effect. B _j ^z0 (t _i ) in the first term on the right side of equation (9) is determined by the value of the spin other than site j at time t _i−1 based on equation (8). B _j ^z (t _i-1 ) in the second term on the right side of equation (9) is s _j ^z (t _i-1 ) / s _j ^x (t _i-1 ) = r _b B _j ^z (t _i-1 ) / B _j ^x (t _i-1 ) and depends on the spin value of the site j (own site) at time t _i-1 . Therefore, the effective magnetic field determination method based on equation (9) is dependent on both the interaction partner and itself at t = t _i-1 and is quantum mechanical in that it depends on both. If the specific effect is seen in the above two-spin example, B _j ^z0 (t _i ) and B _j ^z0 (t _i−1 ), which are opposite signs, are canceled well and the unnatural vibration is eliminated.

In quantum mechanics, the effective magnetic field is determined based on Equation (4). The eigenvalue of σ ^ _k ^z is ± 1. However, in the local field response method, the spin variable s _k ^z is operated so as to take the expected value <σ ^ _k ^z >, so | s _k ^z | ≦ 1. Therefore, in general, the term of Σ _{k (≠ j)} J _kj s _k ^z is underestimated compared to g _j . This situation is not changed by the handling in the expression (9) of the third embodiment. If the Σ _{k (≠ j)} J _kj s _k ^z term is underestimated, the accuracy of the solution decreases. Therefore, the value of s _k ^z will be to normalize the value of g _j to the reference. g _j the factor _{_{_{c i = (Σ k s k}}} z (t i-1) 2 / N) 1/2 and multiplying _{^{_{g j norm (t i) =}}} c i g if _{_j} g _j ^norm (t _i ) And Σ _{k (≠ j)} J _kj s _k ^z terms are almost equal, and the accuracy of the solution is improved. However, let c ₁ = 1 / m, where m is the number of time-axis divisions handled in a discrete manner (accordingly, t _m = τ). This is based on c _i = (Σ _k s _k ^z (t _i-1 ) ² / N) ^1/2 to deal with c ₁ = 0 due to s _k ^z (t ₀ ) = 0. It is.

FIG. 6 shows a flowchart including the above handling. The difference from FIG. 5 is that the procedure 102b is changed to the procedure 102c. In the procedure 102c, handling of the factor c _n = (Σ _k s _k ^z (t _i−1 ) ² / N) ^1/2 is added.

In an actual quantum mechanical spin system, spins always influence each other. That is, the spin σ ^ _j ^z at one site _j affects the spin σ ^ _k ^z at another site k, and conversely σ ^ _k ^z affects σ ^ _j ^z . Therefore, the spin σ ^ _j ^z affects itself via the spin σ ^ _k ^z at site k. This is why the state of a spin in quantum mechanics depends not only on the state of the other party's spin but also on its own spin. The magnitude of the influence on the self through the interaction is proportional to Σ _{k (≠ j)} J _kj ² . By the way, the response function r _b (t) ≠ 1 is due to the interaction. If there is no interaction, r _b (t) = 1. Therefore, δr _b (t) ≡1−r _b (t) ∝Σ _{k (≠ j)} J _kj ² is considered. This can be dependent response quantum effect was average function r _b ⁰ (t) to individual problems. Specifically, 1−r _b ⁰ _mod (t) = (1−r _b ⁰ (t)) Σ _{k (≠ j)} J _kj ² / Σ _{k (≠ j)} ave (J _kj ² ) The response function r _b ⁰ _mod (t) is used instead of the average response function r _b ⁰ (t). Here, ave (J _kj ² ) is an average of J _kj ^{2 of} the problem used in determining r _b ⁰ (t). This improvement allows the response function to be optimized for the actual problem and improves the accuracy of the solution.

Fig. 7 shows a flowchart including the above handling. The difference from FIG. 6 is that the procedure 103 is changed to the procedure 103c.

It can be theoretically shown that the magnitude of the influence on itself through the interaction is proportional to Σ _{k (≠ j)} J _kj ² using third-order perturbation theory. In the third-order perturbation theory, the expected values of effective magnetic field and spin are as follows using the time parameter η≡t / τ.

The subscript (0-3) indicates the inclusion of a perturbation term of 0-3 order. The first term to the third term on the right side of both formulas are the same, but the fourth term is added to formula (14). The fourth term is obtained by multiplying -3 by k = j in the third term, and this term affects the own site j through the site i. It is important that the sign of this term is negative, which partially cancels the first term. Since the second and third terms are zero on average, the average behavior is r _b <1.

The fourth term of equation (14) is proportional to Σ _{i (≠ j)} J _ij ² . That is, 1−r _b (t) ∝Σ _{i (≠ j)} J _ij ² .

Although the method for improving the accuracy of the solution through the quantum mechanical consideration has been described up to the fifth embodiment, the method does not necessarily reach the correct solution. Therefore, in this embodiment, a new auxiliary means is introduced. In the method of this embodiment, the calculation is performed by reducing the transverse magnetic field strength. In this embodiment, this process is repeated n _r times.

FIG. 8 shows an example of an algorithm for that purpose. Procedure 10 is one loop for applying a decrease in the transverse magnetic field. Within each loop, steps 102c and 103c are performed, and i → i + 1 is repeated. if the t _i = t _thn exits the loop to re-set the initial value of the spin proceed to the next loop. The initial value of the spin is s _j ^z (t ₀ ) = 0 in the first loop (n = 1). The second and subsequent loops are determined with reference to the spin values obtained so far.

For example, the initial value of the second loop is set to s _j ^z (t _i ) = s _{j, th1} ^z = s _j ^zfd / m using the optimal solution s _j ^zfd at the end of the first loop (step 111). . For the initial value of the third loop, for example, the optimal solution s _j ^zfd at the end of the second loop is used and the sign of the spin is inverted. That is, s _j ^z (t _i ) = s _{j, th2} ^z = −s _j ^zfd / m. As the initial value of the fourth loop, for example, the value at the end of the third loop is binarized and used. That is, s _j ^z (t _i ) = s _{j, th3} ^z = s _j ^zfd (t _i−1 ) / m. The initial value of the fifth loop is, for example, s _j ^z (t _i ) = s _{j, th4} ^z = s _j using the optimal solution s _j ^zfd at the end of the fourth loop as in the case of the second loop. and ^zfd / m. In this way, if the initial value is changed by various methods, the search space for the solution is expanded and the probability of reaching the correct solution is improved. If the calculation is finished in the fifth loop, n _r = 5.

As shown in step 102c, the strength of the effective magnetic field is

It is. One loop of time range rather than tau t _thn reduced application of transverse magnetic field to implement n _r times - becomes t _{th (n-1).} Coefficients a _n adjusts the amount of change in the lateral magnetic field strength by a factor introduced for that purpose. If a _n (t _thn -t _{th (n-1)} ) = τ, 100% transverse magnetic field is applied. However, it is not always necessary to swing 100%. For example, if a _n (t _thn −t _{th (n−1)} ) = 0.6τ is set, the amount of change in the transverse magnetic field strength can be suppressed to 60%.

In the above example, spin inversion was performed in the third loop. This is because the spin state falling into the local optimum solution is set to the farthest arrangement in the spin arrangement space, so that the global optimum solution can be easily found. In the above example, spin inversion was performed as the initial value of the third loop, but it is also effective to perform spin inversion in the middle of the loop. In this case, time t _i = _tinv may be defined and the spin inversion timing may be controlled.

In the fourth embodiment, g _j (t _i ) is changed to g _j ^norm (t _i ) = c _i in order to contribute the term of g _j (t _i ) and the term of Σ _{k (≠ j)} J _kj s _k ^z equally. g _j was replaced. This process is optimal when the initial spin value is s _j ^z (t ₀ ) = 0. However, in the case of the multi-loop configuration in the present embodiment, the initial value is determined by using the optimum solution s _j ^zfd up to the second and subsequent loops. In this case, since the initial value s _k ^z is an optimal solution candidate, the degree of underestimation for Σ _{k (≠ j)} J _kj s _k ^z is reduced as compared with the case where s _j ^z (t ₀ ) = 0. Therefore, in the second or subsequent loop g _j ^norm (t _i) and multiplying the factor f (≧ 1) as fg _j ^norm (t _i), is effective to expand the contribution of g _j ^norm (t _i) is there. The value of f is determined empirically from about 2 to 50.

This embodiment is shown as an algorithm and can be operated as software on a normal computer or on dedicated hardware. The features of this embodiment are that the operation is relatively simple and the parallelism is high. Therefore, when constructing dedicated hardware, a highly parallel configuration is preferable. If existing hardware is used, it is effective to use highly parallel hardware such as GPGPU (General-Purpose computing (on Graphics) Processing Processing Units). In this embodiment, a configuration example of an apparatus for effectively operating the present invention will be shown.

FIG. 9 shows an example of the computer configuration of this embodiment. FIG. 9 is similar to the configuration of a normal computer, but includes a local field response calculation device 600. The local field response calculation device 600 is a part that specializes in the calculation described in the embodiment 1-6, and other general calculations are performed by the general calculation device 502.

The above configuration may be configured by a single computer, or any part of the main storage device 501, the general arithmetic device 502, the control device 503, the auxiliary storage device 504, the input device 505, the output device 506, etc. You may comprise with the other computer connected with the network.

General operations are operated in the same procedure as a normal computer. Data is exchanged between the main storage device 501 that is a storage unit and the general arithmetic unit 502 that is a calculation unit, and the calculation is advanced by repetition of the data. The control tower at that time is a control device 503 as a control unit. A program executed by the general arithmetic device 502 is stored in the main storage device 501 which is a storage unit. If the main storage device 501 has insufficient storage capacity, the auxiliary storage device 504 that is also a storage unit is used. An input device 505 is used to input data and programs, and an output device 506 is used to output results. The input device 505 includes an interface for network connection in addition to a manual input device such as a keyboard. This interface also serves as an output device.

In the local field response calculation of the present embodiment, N spin variables s _j ^z (t) and N effective magnetic field variables B _{eff, j} ^z (t) are alternately repeated as described in the embodiment 1-6. Ask. The local field response calculation device 600 performs this repeated calculation specialized. On the other hand, calculations other than repetitive calculations such as energy calculation at each time described in the second embodiment are performed by the general calculation device 502.

FIG. 10 shows details of the local field response calculation unit. The local field response calculation device 600 includes a dedicated storage device 601,

calculation units

611, 612, and 613, and registers 621, 622, and 623. The

arithmetic units

611, 612, and 613 are parallel arithmetic devices that operate independently from each other, and each calculates an effective magnetic field B _{eff, j} ^z (t) and a spin variable s _j ^z (t) for each site. The calculation result s _j ^z (t) is stored in the dedicated storage device 601. Information J _ij , g _j , s ₁ ^z (t), s ₂ ^z (t),..., S _N ^z (t) required by the

arithmetic units

611, 612, and 613 are read from the registers 621, 622, and 623. The registers 621, 622, and 623 are dedicated storage units directly connected to the

calculation units

611, 612, and 613, respectively. Thanks to these registers, the

calculation units

611, 612, and 613 can independently perform high-speed processing. In FIG. 10, each register is classified into a region and b region for convenience. J _ij and g _j whose values do not change through the operation are stored in the a region, and the values change with the operation s ₁ ^z (t), s ₂ ^z (t), ..., s _N ^z (t ) Is stored in area b. J _ij and g _j stored in the a area are transferred from the main storage device 501 to the dedicated storage device 601 in advance and further transferred to the

registers

621a, 622a, and 623a. s ₁ ^z (t), s ₂ ^z (t),..., s _N ^z (t) stored in the 621b, 622b, and 623b in the b area are transferred from the dedicated storage device 601. This data transfer is fast because it only sends s ₁ ^z (t), s ₂ ^z (t),..., S _N ^z (t) together, and there is no need for random access. As described above, parallelism and high speed are realized by the configuration of the registers 621, 622, and 623 and the

arithmetic units

611, 612, and 613.

As described in Example 2, the calculation not only repeatedly obtains the effective magnetic field B _{eff, j} ^z (t) and the spin variable s _j ^z (t), but also binarizes s _j ^z (t) at each time, Calculate energy based on it. In this calculation, data is transferred from the dedicated storage device 601 to the main storage device 501, and the general arithmetic device 502 is used. That is, the general arithmetic unit 502 is used except for the repeated calculation of s _j ^z (t) and B _{eff, j} ^z (t). As a result, iterative calculations that require the most calculation time become efficient.

問題 Problems with high difficulty among combinatorial optimization problems belong to NP difficulty. In addition, problems classified as P and problems classified as NP can all be reduced to NP difficult problems. Therefore, almost all combinatorial optimization problems can be solved by solving NP-hard combinatorial optimization problems. The ground state search problem of Equation (1) can also deal with the NP difficulty problem. In this embodiment, the state of the response is shown by taking the maximum cut problem, which is a typical NP difficulty problem, as an example.

The maximum cut problem is a graph theory problem. In graph theory, graph G is composed of vertex set V and edge set E and is written as G = (V, E). The edge e is written as e = {i, j} using two vertices. A graph that is defined by including an orientation in edge e is called a directed graph, and a graph that does not contain a definition of orientation is called an undirected graph. A weight is also defined for the edge e, which is written as w _ij and w _ji . For an undirected graph, w _ij = w _ji . The MAX-CUT problem is a problem of finding a division method that maximizes the sum of the weights of the edges to be cut in the problem of dividing the vertices of the weighted undirected graph G = (V, E) into two groups. If the two undirected graphs after division are G ₁ = (V ₁ , E ₁ ) and G ₂ = (V ₂ , E ₂ ), then the MAX-CUT problem is

It is a problem that maximizes. S _i = 1 with respect to the vertex i∈V _1, if s _j = -1 with respect to the vertex J∈V ₂

Can be written. Since the first term on the rightmost side is a constant if the graph G is determined, the MAX-CUT problem becomes a problem of minimizing Σ _{i> j} s _i s _j . Ising spin glass Hamiltonian

Given in. Therefore, the MAX-CUT problem is equivalent to the ground state search problem of Equation (1) with J _ij = -w _ij and g _j = 0.

The above embodiments are operated on a classic machine, and do not need to be cryogenic and do not need to consider quantum coherence. As a result, a wide range of resources can be used, and electrical circuits can be used. Furthermore, the accuracy of the solution is improved by making the response function dependent on two consecutive times. With these properties, a practical computer that can solve difficult problems with high accuracy is realized.

The present invention is not limited to the above-described embodiment, and includes various modifications. For example, a part of the configuration of one embodiment can be replaced with the configuration of another embodiment, and the configuration of another embodiment can be added to the configuration of one embodiment. In addition, a part of the configuration of each embodiment can be added to or replaced with the configuration of another embodiment, or can be deleted from the configuration of another embodiment.

It can be applied to various data analysis processing.

10 to 300 represents each procedure in the flowchart
501 Main memory
502 General arithmetic unit
503 controller
504 Auxiliary storage device
505 input device
506 output device
600 Local field response calculation device
601 Dedicated storage device
611, 612, 613 arithmetic unit
621a, 621b, 622a, 622b, 623a, 623b registers

Claims

A computer comprising a calculation unit, a storage unit, a control unit, and performing calculation while exchanging data between the storage unit and the calculation unit under the control of the control unit,
N variables s j z (j = 1, 2, ..., N) take a range of -1 ≤ s j z ≤ 1, and a coefficient g j representing a local term and a coefficient J kj ( k, j = 1, 2,…, N)
In the calculation unit, the time is divided into m and discretely calculated from t = t 0 (t 0 = 0) to t m (t m = τ),
At each time t i (i = 1, 2, .., N), variables B j z0 (t i ), B j z (t i ), and s j z (t i ) are determined in order, and time t The initial values of 0 are B j z (t 0 ) = 0 and s j z (t 0 ) = 0, and the above B j z0 (t i ) is B j z0 (t i ) = (Σ k (≠ j) defined by J kj s k z (t i -1) + g j), wherein B j z (t i) by using a parameter u that satisfies 0 ≦ u ≦ 1 B j z (t i) = (1- u) B j z0 (t i ) + uB j z (t i−1 ), B j z (t i ) is multiplied by factor t i / τ and B eff, j z (t i ) = B j z (t i ) · t i / τ
The s j z (t i ) is defined by s j z (t i ) = f (B eff, j z (t i ), t i ) using the function f, and the function f is s j z ( t i ) is defined to be −1 ≦ s j z (t i ) ≦ 1,
Close time step to -1 or 1 the variable s j z as proceeding from t = t 0 to t = t m, finally s j z <0 if s j zfd = -1, s j z> 0 If so, a computer characterized by determining the solution as s j zfd = 1.
For the function f
Using a constant γ , let B eff, j x (t i ) = γ (1 − t i / τ) and tan θ = B eff, j z (t i ) / B eff, j x (t i ) S j z (t i ) is defined by s j z (t i ) = sin θ, and thus the function f is f (B eff, j z (t i ), t i ) = sin (arctan The computer according to claim 1 , wherein (B eff, j z (t k ) / B eff, j x (t k ))}.
Add correction parameters r s and r b for the function f,
Defines θ tanθ = r b · B eff , j z (t i) / B eff, the j x (t i), s j z (t i) = the by r s · sinθ s j z ( t i ) So that the function f is f (B eff, j z (t i ), t i ) = r s · sin {arctan (r b · B eff, j z (t k ) / B eff, 3. The computer according to claim 2, wherein j x (t k ))}.
If s j z (t i ) <0 at each time t i , s j zd (t i ) = −1, and if s j z (t i )> 0, s j zd (t i ) = 1, s If j z (t i ) = 0 then s j zd (t i ) = 0 and H p (t k ) =-Σ k> j J kj s k zd (t i ) s j zd (t i ) -Σ j g j s j zd (t i ) is calculated at each time t i , and s j zfd = s j zd (t i ′ ) at time t i ′ at which H p (t i ) becomes the minimum value The computer according to claim 1, wherein the computer is a solution.
Wherein g j with g j norm (t i) = (Σ k (s k z (t i-1)) 2 / N) defines a 1/2 · g j, the B j z0 (t i) The computer according to claim 1, wherein B j z0 (t i ) = (Σ k (≠ j) J kj s k (t i−1 ) + g j norm (t i )).
Wherein the correction parameter r b defines the δr b ≡1-r b respect, δr b (t) αΣ k (≠ j) J kj 2 and the computer according to claim 3, characterized in that.
Define n r +1 times t thn (n = 0, 1, 2,…, n r ) and n r coefficients a n (n = 1, 2,…, n r ), and t th0 = 0 and t thnr = and tau, wherein B eff, j z (t i ) = B j z (t i) · t i / τ a t th (n-1) ≦ t <B eff in the time domain of t thn, 2. The computer according to claim 1, wherein j z (t i ) = B j z (t i ) a n (t i −t thn ) / τ.
If s j z (t i ) <0 at each time t i , s j zd (t i ) = −1, and if s j z (t i )> 0, s j zd (t i ) = 1, s If j z (t i ) = 0 then s j zd (t i ) = 0 and H p (t k ) =-Σ k> j J kj s k zd (t i ) s j zd (t i ) -Σ j g j s j zd (t i ) is calculated at each time t i , and if H p (t i ) <H p (t i−1 ), then s j zfd = s j zd (t i ) and t <according to claim 7, wherein the determining the spin values s j z (t thn) at the time t = t thn using spin values s j zfd obtained in the time domain t thn calculator.
n q pieces of time t invn (n = 1, 2 , ..., n q) defining a time t i = t invn in B j z0 (t i) = (-Σ k (≠ j) J kj s k ( 2. T i−1 ) + g j ) and B j z (t i ) = (1−u) B j z0 (t i ) −uB j z (t i−1 ) Listed calculator.
The correction parameter r to s and r b, the computer according to claim 3, characterized in that dependent on the t i and the B j z (t i).
An input device, an output device, a storage device, a general arithmetic device, and a local field response arithmetic device;
N variables s j z (j = 1, 2, ..., N) take a range of -1 ≤ s j z ≤ 1, and a coefficient g j representing a local term and a coefficient J kj ( k, j = 1, 2,…, N)
In the local field response calculation device,
The time is divided into m and discretely calculated from t = t 0 (t 0 = 0) to t m (t m = τ)
At each time t i (i = 1, 2, .., N), variables B j z0 (t i ), B j z (t i ), and s j z (t i ) are determined in order, and time t The initial values of 0 are B j z (t 0 ) = 0 and s j z (t 0 ) = 0, and the above B j z0 (t i ) is B j z0 (t i ) = (Σ k (≠ j) J kj s k z (t i−1 ) + g j ), and B j z (t i ) is expressed as B j z (t i ) = (1- u) B j z0 (t i ) + uB j z (t i−1 ), B j z (t i ) is multiplied by factor t i / τ and B eff, j z (t i ) = B j z (t i ) · t i / τ
The s j z (t i ) is defined by s j z (t i ) = f (B eff, j z (t i ), t i ) using the function f, and the function f is s j z ( t i ) is defined to be −1 ≦ s j z (t i ) ≦ 1,
As the time step is advanced from t = t 0 to t = t m , the variable s j z approaches −1 or 1,
In the general arithmetic unit,
The time t i s in each j z (t i) <0 if s j zd (t i) = -1, s j z (t i)> 0 if s j zd (t i) = 1, s If j z (t i ) = 0 then s j zd (t i ) = 0 and H p (t k ) =-Σ k> j J kj s k zd (t i ) s j zd (t i ) -Σ j g j s j zd (t i ) is calculated at each time t i ,
A computer characterized in that s j zfd = s j zd (t i ′ ) at time t i ′ at which H p (t i ) becomes the minimum value is set as a final solution.

calculator.
A calculation method using a computer that includes a calculation unit, a storage unit, a control unit, and performs calculation while exchanging data between the storage unit and the calculation unit under the control of the control unit,
N variables s j z (j = 1, 2, ..., N) take a range of -1 ≤ s j z ≤ 1, and a coefficient g j representing a local term and a coefficient J kj ( k, j = 1, 2,…, N)
In the calculation unit, the time is divided into m and discretely calculated from t = t 0 (t 0 = 0) to t m (t m = τ),
At each time t i (i = 1, 2, .., N), variables B j z0 (t i ), B j z (t i ), and s j z (t i ) are determined in order, and time t The initial values of 0 are B j z (t 0 ) = 0 and s j z (t 0 ) = 0, and the above B j z0 (t i ) is B j z0 (t i ) = (Σ k (≠ j) J kj s k z (t i−1 ) + g j ), and B j z (t i ) is expressed as B j z (t i ) = (1- u) B j z0 (t i ) + uB j z (t i−1 ), B j z (t i ) is multiplied by factor t i / τ and B eff, j z (t i ) = B j z (t i ) · t i / τ
The s j z (t i ) is defined by s j z (t i ) = f (B eff, j z (t i ), t i ) using the function f, and the function f is s j z ( t i ) is defined to be −1 ≦ s j z (t i ) ≦ 1,
Close time step to -1 or 1 the variable s j z as proceeding from t = t 0 to t = t m, finally s j z <0 if s j zfd = -1, s j z> 0 If so, a calculation method characterized by determining the solution as s j zfd = 1.
Wherein g j with g j norm (t i) = (Σ k (s k z (t i-1)) 2 / N) defines a 1/2 · g j, the B j z0 (t i) 13. The calculation method according to claim 12, wherein B j z0 (t i ) = (Σ k (≠ j) J kj s k (t i−1 ) + g j norm (t i )).
Define n r +1 times t thn (n = 0, 1, 2,…, n r ) and n r coefficients a n (n = 1, 2,…, n r ), and t th0 = 0 and t thnr = and tau, wherein B eff, j z (t i ) = B j z (t i) · t i / τ a t th (n-1) ≦ t <B eff in the time domain of t thn, 13. The calculation method according to claim 12, wherein j z (t i ) = B j z (t i ) a n (t i −t thn ) / τ.
For the function f
Using a constant γ , let B eff, j x (t i ) = γ (1 − t i / τ) and tan θ = B eff, j z (t i ) / B eff, j x (t i ) S j z (t i ) is defined by s j z (t i ) = sin θ, and thus the function f is f (B eff, j z (t i ), t i ) = sin (arctan (B eff, j z (t k ) / B eff, j x (t k ))}
Add correction parameters r s and r b for the function f,
Defines θ tanθ = r b · B eff , j z (t i) / B eff, the j x (t i), s j z (t i) = the by r s · sinθ s j z ( t i ) So that the function f is f (B eff, j z (t i ), t i ) = r s · sin {arctan (r b · B eff, j z (t k ) / B eff, j x (t k ))}
Wherein the correction parameter r b defines the δr b ≡1-r b respect, δr b (t) αΣ k (≠ j) J kj 2 and calculation method of claim 12, characterized in that the.