Fully and Partly Divergence and Rotation Free Interpolation of Magnetic Fields

A new interpolation method for rotation and divergence free fields is presented. It is based on a suitable choice of a tricubic interpolation scheme and reaches an accuracy of third order in grid size (Δx). With the interpolation method it is possible to increase the accuracy with a factor of grid size/distance with respect to the trilinear interpolation method using exactly the same data points. Simulations for several distances of dipoles (r) to the interpolation area show that the maximum relative deviation is approximately 3(Δx/r) ppm.


Introduction
For many applications it is needed that a vector field is accurately known with sufficient spatial resolution.For instance for the calculation of neutron spin precession [1], in magnetic resonance imaging [2] or in space [3].In case these magnetic fields can be calculated, it is possible that these calculations are very time consuming and it might be impossible to do these calculations with the required spatial resolution.In case these magnetic fields have to be measured, it is also possible (and even likely) that these measurements are very time consuming and the required spatial resolution can not be obtained.In these cases one must resort to interpolation methods.It is important to have an interpolation method that is accurately enough, even when the grid size of the calculated or measured data is large with respect to the required resolution.This calls for higher order interpolation schemes.The use of a divergence and rotation free interpolation can be advantageous in two respects.First, by using these properties of the field it is possible to increase the accuracy of the interpolation.Second, in some cases if the interpolated field is not rotation and divergence free, this leads to erroneous results which have no physical meaning.This is the case in for instance magneto-hydrodynamics [4].An important example of a rotation and divergence free vector field is the magnetic field in a conduction current free region in a homogeneous material.The magnetic field is always divergence free and when there are no conduction currents, the field is also rotation free.In literature various methods exist based on a trilinear [5], triquadratic, or even tricubic [6] interpolation schemes.The disadvantages of these methods are either their limited accuracy (trilinear and triquadratic) or the large number of grid points needed (triquadratic and tricubic).Further, in general, these interpolation methods result in a field inside the grid that has a finite divergence and rotation, possibly resulting in unphysical results.Here, a method is presented that gives a rotation and divergence free interpolation of a vector field inside a 3D-rectangular grid, with an accuracy that is third order in grid size.

Method
Lets consider a rectangular box of homogeneous material (with linear constitutive equations) with sides Δx, Δy and Δz (see also Figure 1).The Maxwell equations require that the magnetic induction is always divergence free.In case of the quasi-static limit and without conduction currents the magnetic induction should also be rotation free.As the material is homogeneous and linear, these conditions also hold for the magnetic field, H .
A rotation free field can always be described as the gradient of a scalar quantity, so let us define a scalar field, G so that because of the fact that 0   H  .It is assumed that G is a mixed polynomial of x, y and z , , , , where i, j, k are elements of   0,1, 2, , 1 m   .In the following it is explained which choices where made to find the coefficients , , i j k .0,0,0 is taken 0 as only the gradient of G is of importance.Condition (2) yields where i', j', k' are elements of .Hence there are equations (note that i', j' and k' have m − 2 values and i, j and k have m values).This can be rewritten as where all coefficients with any index i, j or k larger than m − 1 should be taken 0. Further, all coefficients with any two indices larger than m − 3 should also be 0 to get a divergence free field.These are 7 equations for m = 3 and 12(m − 4) + 8 equations for m > 3.
For m = 3 there are 12 independent coefficients.It seems a logical choice to use 4 points of the rectangular grid.The complete rectangular box can be covered by dividing it in 6 tetrahedrons, A to F. Each point in the box lies in exactly one of these tetrahedrons (see Figure 1).On each edge of the tetrahedron the three components of the magnetic field are known, hence there are 12 independent values that determine the coefficients of G.All regions should contain a body diagonal to get independent values of locations of the magnetic field components.In any other case no solution can be found, because the matrix as defined in Equation ( 4) will be singular, as some rows of the matrix will be identical.It is not directly obvious, but it can be derived by trying all possibilities.

Q
The scalar G is a polynomial with a maximum of second order in x (or y, z), so that the magnetic field component in the x (or y, z) direction is only a function of first order in x (or y, z).The interpolation will only be accurate up to first order.This is a similar accuracy as is obtained with simple trilinear interpolation of each separate component [5].
For m > 3, there are (m − 2) 3 + 6 (m − 2) 2 + 12 (m − 4) + 8 equations and the number of coefficients is m 3 − 1, so that the number of free coefficients is 23 and independent of m.For m = 4 the scalar G is a polynomial with a maximum of third order in x (or y, z), so that the magnetic field is of second order in x (or y, z) so that the interpolation can be accurate up to second order.For m = 4, G can be written as where the elements of A and X are given by the first 23 numbers in Table 1.The magnetic field is given by Equation ( 1), so that and similar for the y and z components.As A contains 23 coefficients, 23 magnetic fields components are needed to be able to determine these coefficients.Of course as only 23 components are needed one is not used.
As soon as the interpolation gets close to this point it will get less accurate as one component is inferred from the other components.This can be prevented by dividing the box in two regions similar to Figure 2 so that the interpolation will always be at sufficient distance from this point.
Another option is to add one additional coefficient.The disadvantage is that in such a case the exactness of the divergence free interpolation is given up.It has been found that this method gives slightly more accurate results for the simulations considered here.Now, A and X contain 24 elements that are given by the complete Table 1.
When is a vector with 24 elements, determined by the values of the field components at the corners of the rectangular box, then it is possible to write where is a matrix connecting the elements of both vectors.Each element of a row of is determined by the corresponding element of Equation ( 3) where i is to be taken over the 8 corners of the rectangular box, corresponding with coordinates (x i , y i , z i ).Equa-tion (3) can be inverted to find the elements of and Equation ( 3) can be used to perform the interpolation.
For a regular grid matrix is fixed and has to be determined only once.One can even perform the inversion of the matrix analytically, to speed up the process.For a cubic grid and Δx = 1, the matrix reduces to a very simple one and is shown in Table 2.The advantage of this method over the one presented by Lekien [6] is that no derivatives need to be determined at the corners and that matrix has only 24 components instead of 96.

Simulation and Results
The accuracy of the interpolation scheme can be tested by means of the magnetic field of a dipole at a certain distance from the rectangular box that is considered.The calculated magnetic field at a certain point in the box can be compared to the one that is interpolated using the magnetic field points at the corners of the box.This has been done for above interpolation schemes (m = 3, m = 4 and 24-element) for a cubic grid and magnetic dipoles directed in the x, y or z direction and positioned at the x-axis at a certain distance.The magnetic field of such a dipole was calculated according to where µ 0 is the magnetic permeability of vacuum, is the vector from the dipole to the interpolation point and is the magnetic dipole vector.For a magnetic dipole in the x-direction only the x-component of is non-zero (and similar for the y-and z-directions).Examples of the results for the 24-elements interpolation scheme are shown in  the magnitude of the calculated magnetic field as function of x and y for z = Δx/2.The distance between the box and the magnetic dipole is 10Δx (left side) or 100Δx (right side) and the results are shown for three different dipoles (from top to bottom directed in the x, y and z direction).The relative deviation is already quite small for a dipole at a distance of 10 times the grid size and reduces another factor of 1000 when the distance is increased by a factor of 10.The graphs show that the relative deviation is not constant but depend on the properties of the field.The interpolation coefficients are determined by the values at the corners of the cube, but the deviations do not approach zero at the corners.This is because a slice through the center of the cell (as used in the plots) does not contain any corners.This slice was chosen to show the larger deviations.The interpolation is accurate to third order in grid size.Further simulations for several distances of the dipoles, r to the interpolation area show (see Figure 4) that the maximum relative deviation is approximately 3(Δx/r) 3 ppm.If the interpolation scheme for m = 4 is used with 23 elements, the maximum relative deviation that occurs in this case is approximately 5(Δx/r) 3 ppm.The interpolation scheme for m = 3 or the trilinear interpolation have maximum deviations of approximately 4.5(Δx/r) 2 ppm and 1.5(Δx/r) 2 ppm respectively.Note that the latter interpolation schemes are only accurate in second order of the grid size.Hence, by using the properties of the magnetic field (divergence and rotation free) it is possible to increase the accuracy of the interpolation with a factor of Δx/r with respect to the trilinear interpolation using exactly the same data points.

Conclusion
A new interpolation method uses the properties of the magnetic field in homogeneous linear materials without conducting current (i.e.rotation and divergence free).The method needs 24 magnetic field components, which are provided by the 24 components available at the 8 corners of a rectangular grid.The interpolation accuracy obtained is third order in grid size.With the interpolation method it is possible to increase the accuracy with a factor of grid size/distance with respect to the trilinear interpolation using exactly the same data points.The interpolation method can be checked for other field shapes as used here to check the validity of the results further.
The number of free parameters for a divergence and rotation free field derived from a potential represented by a polynomial with m > 3 is always 23.It would be inter-esting to create a similar interpolation scheme for m > 4, to find out the influence on the accuracy.

Figure 1 .
Figure 1.Grid points with different regions A to F. Note that the different regions rotate around the body diagonal.

Figure 3 .Figure 2 . 2 + 2 .Figure 3 .
Figure 2. Different regions A and B in 23 coefficients interpolation.In region A the z-component of the magnetic field B 111 is ignored and in region B the z-component of B 000 .The dividing plane is 2zΔz + 2yΔy + 2xΔx = Δx 2 +Δy 2 + Δz 2 .

Figure 4 .
Figure 4. Relative deviation of interpolation (24 elements) as function of distance to dipole for dipoles in x, y and z direction.The full black line is the maximum limit given by 3(Δx/r) 3 ppm, r is the distance to the dipole.The (upper) red line indicates the limit for interpolation with m = 3, the long-dashed-green line the limit for trilinear interpolation and the short-dashedblue line the limit for interpolation with m = 4.