Magnetic Signature Identification using Pareto optimization

Christopher Riley, Engineering Manager, Cobham Technical Services

Abstract

Determining the topology of a naval vessel from its measured magnetic signature is a classic inverse problem – what is the source of the observed response? This paper explores an alternative to classical inverse methods with a technique that is used in design optimization. Pareto optimization is commonly used to synthesise the design of new equipment, where there is a requirement to explore a design space to reach multiple objectives. In the inverse application, the object has already been designed and its measured performance is known. However, its exact geometry or, even, topology and its location are unknown. The paper shows that Pareto optimization can be successfully used to accurately determine overall dimensions of a simplified model of a submarine pressure hull from its magnetic signature. It also shows that Pareto optimization applied to a simple model consisting of two hollow triangular prisms can be used to give reasonable knowledge of the overall dimensions of either a surface ship or a submarine.

1. Introduction

Magnetic signatures of naval vessels occur because of the interaction of ferromagnetic and conducting materials in the platform with the Earth’s magnetic field [1], [2]. There are various types of signatures that can occur. The permanent signature, caused by realignment of domain walls beyond their reversible energy, occurs when the vessel maintains a constant heading in Earth’s field over a considerable period of time, is subjected to a violent shock or undergoes very large changes in mechanical stress. The induced signature results from the instantaneous response of the magnetization due to the heading at a particular time. Permanent and induced signature are often of similar magnitude, although most vessels undergo a de-perming process to try to reduce the former. A smaller signature, called the eddy current signature, also occurs due to the motion of the vessel – especially its side-to-side rocking motion [3]. The motion of the ferromagnetic and other conducting materials in the Earth’s field induces eddy currents and these, in turn, cause a reaction field. In addition, all vessels also include equipment for cathodic protection, which creates a DC electric field surrounding the platform and, in the sea, causes currents to flow and return through the hull. These currents produce the so-called corrosion-related magnetic signature field.

The consequence of this is that sea based mines have been developed that recognise the magnetic anomaly caused by the signature and detonate. Equally, there is considerable activity in the measurement of magnetic anomalies – both from detection sub-sea and airborne – in an effort to establish which platform is in the vicinity. To overcome this, countermeasures, such as de-perming and degaussing, are often used to reduce signature. Nevertheless, magnetic anomaly detection (MAD) is still widely used and this has led to much interest in the solution of inverse problems to determine the platform causing the anomaly. Broadly stated, an inverse problem can be described as determination of an unknown source from a known response and there are many classic methods that can be applied [4].

In this paper, a method normally used in design synthesis [5], [6] has been applied to the MAD inverse problem. Pareto methods are commonly used when solving design optimization problems with multiple, usually competing, objectives. For example, produce the maximum level of central field maintaining the best possible homogeneity over a defined volume using less than a pre-defined maximum current limit in a fixed number of coils. Continuous design variables, such as current in and position of each coil, result in a multi-dimensional design space to search, subject to the constraint on the current limit. The two objectives – homogeneity over the largest volume and maximum field – are likely to be in conflict. Homogeneity will tend to spread the position of the coils and make the current in each similar, whereas maximum central field will tend to concentrate the coils towards the centre with variable currents in each coil. One method of solving this multi-objective problem is to assign a weight to each objective and add all the weighted objectives to form a single objective – leading to only one optimal design. In contrast, Pareto methods produce a multidimensional optimization space with a family of optimal designs that are on the Pareto optimal “front”. These are designs for which the Pareto algorithm can find no other designs for which all of the objectives are improved. A simple two objective result is shown in figure 1. In this case, minimization of both objective functions leads to a family of six Pareto optimal designs. Some are better at minimizing one objective, others the second objective, while some offer different levels of compromise. The infeasible solutions also shown on the optimization space map are designs that have been investigated but fail one or more of the constraints applied. Feasible solutions satisfy the constraints but are “dominated” by the Pareto optimal designs.

The Pareto method has also been quite successfully used to determine size and position of buried objects [7], such as oil and gas pipelines. In this work, the topological shape of the object (the pipeline) was already known, but its position below ground level relative to an electromagnetic AC coil at the surface and the diameter and wall thickness of the pipeline were unknown. However, since the impedance of the coil was known, the objective functions minimized could be defined as the difference between the computed and measured resistive and inductive components of the impedance. It was shown that designs on the Pareto optimal front tended to “cluster” in design variable space i.e. several designs on the front had similar parameter values. At one particular AC frequency, several of these clusters may form but, when more than one frequency is used, only one cluster appears at all frequencies. It was shown that this unique cluster had a reasonable match to the true dimensions and position, although the thickness of the pipe wall was difficult to determine.

In this paper, the same technique has been extended to measured induced signatures. Field values on a two dimensional grid below the vessel (as would typically be obtained from a range) are used in the objective functions. Two branches of the research are reported here: identification of a known source topology (a hollow steel tube with hemispherical end caps, but unknown dimensions and wall thickness) and an unknown topology (signatures from models of either a surface ship or a submarine). In all cases, the parameterized model that is used in the Pareto optimization is represented by a 2D thin-plate surface in a non-linear magnetostatic simulation using the Opera software [8], [9]. The 2D thin-plate surface simplifies model and mesh construction by eliminating the requirement to handle high-aspect ratio volume structures (such as deck and bulkhead plates). It also gives considerable benefit in reducing simulation requirements (memory and run time) for the same accuracy, compared to the equivalent model constructed with the actual volumes.

2. Inversion study with a known topology

The known topology study was performed using magnetic simulation results from a simple representation of the pressure hull of a submarine, as shown in figure 2. This consists of a hollow steel tube of length 45 metres, a radius of 4 metres and a steel thickness of 150 mm, with hemispherical end caps made from the same thickness steel. Two sets of signature magnetic field results were obtained for this geometry at a depth of 25 metres below the tube axis. In one case, Earth’s field was aligned along X (longitudinal) with a component of -40 A/m and Z (vertical) component of +50 A/m. In the second case, the Z component was the same but the horizontal component of -40 A/m was now applied in Y (athwartships). Figure 2 also shows a typical result showing the Z component signature from the second case. The field is sampled on a 100 x 100 metres XY regular grid at a spacing of 2.5 metres.
Two types of minimization criteria were examined for the optimization at each calculated design instance during the inversion. In one type, the magnitude of the difference between the measured and calculated integral of a particular field component over the 100 x 100 metres grid was used, resulting in two objectives to minimize. (Note that the integral in the 3rd component direction will be zero because of the field and geometry symmetry – so any difference is just numerical noise). In the other type, the difference between the field components at discrete points on the grid was minimized to give 4 objectives. The points chosen were where the maximum positive and negative signature field occurred in one of the applied field component directions for the measured values. In all cases, the range searched for the three design parameters was:

• 25 ≤ length ≤ 60 metres
• 3 ≤ radius ≤ 5 metres
• 100 ≤ thickness of steel ≤ 200 mm

Although this inversion is not dissimilar to the previous work undertaken on buried objects, the behaviour is not entirely similar. For the integral objectives, there were only 4 Pareto optimal designs found when the Earth’s field was applied in XZ and all were in the same cluster very close to the true dimensions for length and radius, as shown in figure 3a, where the parameter values have been normalized against the real values. However, finding the thickness of the steel was more varied, although the mean value of the 4 designs agreed to about 5%, as shown in table 1. As shown in figure 3b, when the Earth’s field was applied in YZ and integral objectives were used, only a single Pareto optimal design was found. This again showed very good agreement in length and radius prediction and was 4% in error for the steel thickness. Consequently, there is no need to identify clusters in each set of results and determine which cluster exists in both.

A more complex picture emerges when the objectives are based on signature field values at discrete points, as shown in figures 4(a) and (b) for XZ and YZ applied field respectively. For XZ applied field, in all Pareto optimal designs the length and radius are still clustered close to the correct value. This can also be seen in table 2 from the ~1% standard deviation of both these values. As with the integral objectives, prediction of the steel thickness is varied, although designs 4, 5, 6 and 9 are all within 10% of the correct value.

For the YZ applied field, as shown in figure 4(b), there is a cluster of designs within 1% of the correct values for both length and radius, comprising six designs: 1, 3, 4, 6, 7 and 11. There are also three clusters each containing 2 similar designs for length and radius: 10 and 13, 9 and 14, 2 and 8. These clusters tend to over predict one of the length or radius parameters and under predict the other one. However, when the field was applied in XZ, as discussed in the previous paragraph, no designs with the values in the small clusters emerged. Consequently, it is still possible to identify which of the clusters in the YZ applied field results is showing the approximate correct values. As with all other optimizations, the prediction of steel thickness is much more varied. However, the average of the six designs in the main cluster is 157 mm – only 4.7% in error.
Tables 1 and 2 also show that the cluster around the correct dimensions tends to dominate the overall results. Although smaller clusters may form, they do not provide sufficient weight, when the average value of all designs is taken, to disturb the average value. In three of the four optimizations, only the cluster in the vicinity of the actual values emerged, while in the Discrete Values YZ Applied Field results, the average length for the six designs in the dominant cluster is 45.47 metres and, for the radius, the average 3.994 metres. These are not significantly better than the mean for all designs (45.92 and 3.995 metres). Consequently, in the subsequent inversion from signatures of an unknown topology, it was decided to use average values from all the Pareto optimal designs rather than identify clusters. It was also decided to only use the integral objectives, as these perform more reliably (as shown from the standard deviation values in table 2).

3. Inversion study with unknown topologies

In the previous section, Pareto optimization was successfully used to determine, with a reasonable degree of accuracy, the length and radius of a structure with a known topology (a hollow tube with end caps). Material thickness is more problematic with a much wider dispersion of results, probably due to the fact that the material is not magnetically saturated and significant variation in magnetization does not occur.

However, when MAD results are obtained, the type of vessel creating the signature is often unknown. Is it a small vessel nearby or a larger vessel at some distance? Is it a surface ship or a submarine? To investigate the possibility of using Pareto optimization as an inversion tool in this situation, it was decided to investigate models that might, at least, be able to predict overall dimensions of an unknown topology. Two unknown topologies were chosen: a surface ship, as shown in figure 5, and a submarine, shown in figure 6. The finite element models for these vessels have been cut away to reveal some detail of their internal structure. As can be seen, the ship contains both internal decks and bulkheads, while the submarine is double hulled. The overall dimensions of these vessels are also significantly different, as shown in table 3.

The number of geometric parameters that would be required to describe all the important dimensions of either of these platforms would produce a very complex, multi-dimensional design space for optimization. Realistically, 10 to 12 dimensional space is probably the limit of what can be undertaken, as will be shown later. Equally, without the constraint on how many dimensions are practically viable, it would be difficult to develop a single parameterized model that could model either topology. Consequently, the models used to predict the unknown topologies only range from 4 to 12 dimensions and use generic shapes: ellipsoids, spheres, prisms and cuboids. In all the studies, the ranges for the dimensions used should be chosen to allow accurate prediction of length, beam and vertical height for either the surface ship or the submarine after a successful inversion, although not all the simplified models use these dimensions explicitly. All studies were also limited to solving a maximum of 250 finite element models for the simplified model and 2 hours total computation time.

a. Hollow Ellipsoid

The ellipsoid is characterized by 4 dimensions: 3 half-axis lengths and the steel thickness. The ranges for the dimensions are shown in table 4.

Table 5 shows the average value for the values of length, beam and vertical height obtained using XZ and YZ applied Earth’s field for both ship and submarine signatures. The values are normalized to the dimensions shown in table 3. Table 5 also shows the number of Pareto optimal designs found, from which the average values were obtained.
It can be seen that prediction with the ellipsoid model does not produce a satisfactory result for the dimension that is orthogonal to the two components of applied field. Both the submarine and ship lengths are not predicted accurately when YZ field is applied, and the beam result for the ship is very poor for the XZ applied field. It was assumed that this results from the symmetry in the signature of the ellipsoid imposing constraints that do not exist in the measured signatures. As a result, in the other simple models reported later, it was decided that the Earth’s field would be applied in a single simulation that would avoid any symmetrical measurements. The horizontal component of -40 A/m was moved to be at 45 degrees to both the X and Y directions, giving a significant non-zero integral in all three components. In general, the prediction for the submarine dimensions is better than for the ship. This probably is a result of the geometry of the submarine being closer to an ellipsoid shape.

b. Hollow spheres

5 hollow steel spheres were defined situated on the X-axis at centres x = -30, -15, 0, 15, 30 metres, each with their radius as a variable in the range 0.5 to 15 metres. Hence, a five dimensional design space is searched to obtain the Pareto optimal designs. This proved quite unsuccessful. For example, figure 7 shows the configuration obtained from the average parameters using the signature field from the ship model.

c. Two hollow equilateral triangular prisms

5 dimensions were used to define the geometry of the two prisms, shown in figure 8. These were:

• The lengths of the two prisms (2 parameters)
• The radii of the circumscribed circles of the equilateral triangles defining the prisms (2 parameters)
• The vertical height of the hull prism centre line (1 parameter)

Table 6 shows the range for each parameter

Table 7 shows the average value for the values of length, beam and vertical height obtained using XYZ applied Earth’s field for both ship and submarine signatures. The X and Y components are equal, such that their combined field is 40 A/m. The values are normalized to the dimensions shown in table 3. The table also shows the number of Pareto optimal designs found, from which the average values were obtained. As can be seen, the overall dimensions are predicted to within 15% of the real dimensions for both the submarine and the ship. Figures 9 and 10 show the predicted compared to the actual geometry for the submarine and the ship respectively, with both the predicted and actual geometry drawn at the same scale.
d. Combination of three hollow prisms and a hollow cuboid
The final simple approximation was more complex than the previous ones and, to some extent, is a closer approximation to the shape of a real vessel, as shown in figure 11. Twelve parameters were used to describe this approximation.
• The vertical height and the radius of the circumscribed circle for each prism (6 parameters)
• The length of each prism (3 parameters)
• The length, width and height of the cuboid (3 parameters)
Table 8 shows the range for each parameter
Table 9 shows the average value for the values of length, beam and vertical height obtained using XYZ applied Earth’s field for both ship and submarine signatures. The X and Y components are equal, such that their combined field is 40 A/m. The values are normalized to the dimensions shown in table 3. The table also shows the number of Pareto optimal designs found, from which the average values were obtained. As can be seen, this simple approximation does not perform as well as the previous two prism model, with the worst case being the normalized beam of the submarine which is more than 100% in error. Figures 12 and 13 compare the actual and predicted geometry for the submarine and ship respectively.

4. Discussion

As previously mentioned, the results reported in section 2, using the known topology example, show that objectives for minimization derived from the integrals of many values of magnetic field give improved inverse results from objectives derived from fields at discrete points. This has also been seen in other applications of Pareto optimization to solve inverse problems in electromagnetics, not reported in public literature. An example is the actual placement of a collection of permanent magnet blocks in a magnetron compared to the engineering design positions, where the integrals of measured and predicted fields at the surface of the magnetron target were used in the objective functions. The behaviour of the magnetron is critically dependent on the magnetic field and tolerance in the assembly of the magnet blocks gave rise to higher erosion of the target in some areas than expected.
The integrals used in the magnetron work were subtly different from those used in the signature work reported here. For the signatures, a typical objective function to be minimized, f_x, was:

formula-1

where S is the area of the grid over which samples have been measured. In the magnetron work, the objective function was defined as

formula-2

The second of these expressions is usually to be preferred, as it will tend to minimize the error between calculated and measured results at every sample point on the grid. However, it has the disadvantage that the measured result at every sample point must be accessed every time the objective function is calculated. In the first expression, the integral of the measured results need only be evaluated prior to the optimization and its value used. However, the first expression is potentially dangerous as the integral of the calculated results may be similar to the measured, while the discrete values have many large positive and negative errors that are tending to cancel. The first expression is probably successful in the signature modelling, as the calculated signature fields from the models produced during the optimization always tend to have the same shape over the grid space as the measured results. Consequently, the possibility of significant cancelling errors at the discrete points is much reduced.
The results presented in section 3 show the fine line that it is necessary to draw between a sufficient number of variables to reasonably capture possible geometry and choosing too many variables. As could be seen from the three prism and cuboid results, in the time / maximum designs limit imposed on all four simple models for effective comparison, it was not possible to adequately sample twelve dimensional space to obtain a statistically significant set of Pareto optimal designs. For the submarine, the magnitude of the standard deviation on the 12 variables ranged from 30 to 80%. For the ship, the range was even larger (20 to 110%), although the majority were around 30 to 50%. With only 6 and 10 Pareto optimal designs for the submarine and ship respectively, a dominant cluster is not emerging.
However, it is also shown that choosing a topology that has some similarity to the expected type is of benefit. While the dimensional space of the ellipsoid and spheres are smaller, they are not sufficiently close to the real topology for good results. In this respect, the two prism model is the most successful with only a five dimensional design space and a topology that has the fundamentals of the ship or submarine (hull and superstructure or sail).
Further work should investigate whether the objective functions used could be improved (as discussed) and whether a better, low dimensional space model than the two prisms can be discovered.

5. Conclusions

The results presented in this paper show that Pareto optimization is an effective method of inversion from a measured MAD signature to determine overall geometry. For a simple known topology, it is possible to obtain very accurate prediction of overall dimensions. Prediction of material thickness is more problematic, since typical Earth’s field values are not sufficiently high to cause magnetic saturation.

For inversion when the topology is complex and unknown, the results have also shown that it is possible to predict overall dimensions of a platform using Pareto optimization. The model chosen should be sufficiently simple in variable dimension space to allow a dominant cluster of results to emerge but have some characteristics that reflect the type of topology for the expected sources. In this respect, a design consisting of two equilateral triangle cross-section hollow steel prisms requiring five variables for its definition was the most successful model examined so far. Prediction of length, beam and vertical height to within 15% of true dimensions, for both a 95 metre long surface ship and a 41 metre long submarine, was obtained. The same ranges for the variables were used in both cases.

The benefit of obtaining measurements when the Earth’s field is in a direction that produces a significant non-zero longitudinal, athwartships and vertical integral of the signature is also shown.
Acknowledgements

The author thanks Mr. Venkatesu Samala, Senior Application Engineer at ICON Design Automation, Bangalore, India for the creation of the submarine Opera finite element model. He also thanks his colleagues at Cobham Technical Services Opera Software for suggestions of simplified model topologies.

References:
[1] John Holmes, “Exploitation of a Ship’s Magnetic Field Signatures – Synthesis Lectures on Computational Electromagnetics”, published by Morgan & Claypole, 2006, ISBN: 9781598290745
[2] John Holmes, “Reduction of a Ship’s Magnetic Field Signatures – Synthesis Lectures on Computational Electromagnetics”, published by Morgan & Claypole, 2008, ISBN: 9781598292480
[3] Paula Y. Zivi, Amy LeDoux, Robert D. Pillsbury, Jr., “Electromagnetic Signature of a Ship Rolling in Earth’s Magnetic Field”, presented at Compumag 2001, Evian, France
[4] P. Neittaanmäki, M Rudnicki, A. Savini, “Inverse Problems and Optimal Design in Electricity and Magnetism”, Oxford Science Publications, 1996, ISBN: 0-19-859383-X
[5] Glenn Hawe, Jan Sykulski, “Scalarizing cost-effective multi-objective optimization algorithms made possible with kriging”, COMPEL Journal, vol. 27, no. 4, 2008
[6] G. I. Hawe, J. K. Sykulski, “An enhanced probability of improvement utility function for locating Pareto-optimal solutions”, Proceedings of COMPUMAG 2007, pp. 965 – 966
[7] C. P. Riley, “Solving Inverse Electromagnetic Sensing Problems using Pareto Optimization”, Sensor Letters, Vol. 11, No. 1, January 2013
[8] C. S. Biddlecombe, C. P. Riley, “Improvements to finite element meshing for magnetic signature simulations”, Proceedings of Marelec 2015, poster 7.
[9] www.OperaFEA.com