Linköping University | Department of Physics, Chemistry and Biology Master of Science Thesis, 30 credits |Applied Physics Spring term 2019 | LITH-IFM-A-EX--19/3635--SE

### Modelling the Magnetic Influence

### of a Jet Aircraft

### A study on the magnetic interference of an

### aircraft configuration and its effect on a

### magnetometer

**Erik Sandlund **

Examiner: Tobias Hansson, IFM, Linköping University Supervisors: Alexander Hutchings, Saab Surveillance

**Datum **

Date

2019-06-20

**Avdelning, institution **

Division, Department

Division of Theoretical Physics

Department of Physics, Chemistry and Biology Linköping University

SE-581 83 Linköping, Sweden

**URL för elektronisk version **

**ISBN **

**ISRN: LITH-IFM-A-EX--19/3635--SE **

_________________________________________________________________
**Serietitel och serienummer ** **ISSN **

Title of series, numbering ______________________________

**Språk **
Language
Svenska/Swedish
Engelska/English
________________
**Rapporttyp **
Report category
Licentiatavhandling
Examensarbete
C-uppsats
D-uppsats
Övrig rapport
_____________

**Titel ** **Modellering av ett jetflygplans magnetiska påverkan **

Title Modelling the Magnetic Influence of a Jet Aircraft

**Författare ** **Erik Sandlund **

Author

**Nyckelord**

Keywords Magnetic fields, magnetic interference, electromagnetic field theory, simulations, engineering, finite-element method, jet engines, aircraft, submarines

**Sammanfattning**

Aircraft have been used for the detection of submarines since World War II. The basic concept is to attach a sensor to the back of an aircraft. Since the aircraft is a moving metallic object, it is bound to generate a great deal of interference. Because of this, mathematical models and software have been developed to help filter out this interference and thus make the detection of the submarine easier. Normally, the engines of the aircraft are placed on the wings, quite far away from the sensor. However, for a maritime patrol system in development, the jet engines are placed at the rear of the airframe, generating the necessity to study whether or not they affect the performance of the sensor, which is the purpose of this thesis.

Several models were created, tested and simulated for the airframe and jet engines. One of each of these were then combined to create a simulation model for the complete aircraft. A jet engine model that included rotating machinery -- a possible source of magnetic interference -- was also created, but could not be added to the model for the complete aircraft. The magnetic interference was mathematically compensated for, removing the static interference, but not the interference during manoeuvres. The jet engine part of the complete aircraft model did not seem to generate a significant amount of magnetic interference compared to the airframe. An electric dipole, representing a submarine, was then added to the simulation. The data from that simulation was put through the mathematical model and distortions of a few~nT were noticeable during straight courses. The jet engine model that included rotating machinery yielded different results compared to the jet engine model in the complete aircraft model. They seemed to contain signals of higher frequency, which were however not detected by a frequency domain study or present during straight courses. It was thus concluded that using this particular engine model the submarine could probably still be detected if the course of the aircraft was kept straight, though further research is needed with more advanced models for the engine, in particular with regards to the rotating machinery.

**Abstract**

Aircraft have been used for the detection of submarines since World War II. The basic concept is to attach a sensor to the back of an aircraft. Since the aircraft is a moving metallic object, it is bound to generate a great deal of interference. Because of this, mathematical models and software have been developed to help filter out this interference and thus make the detection of the submarine easier. Normally, the engines of the aircraft are placed on the wings, quite far away from the sensor. However, for a maritime patrol system in development, the jet engines are placed at the rear of the airframe, generating the necessity to study whether or not they affect the performance of the sensor, which is the purpose of this thesis.

Several models were created, tested and simulated for the airframe and jet en-gines. One of each of these were then combined to create a simulation model for the complete aircraft. A jet engine model that included rotating machinery – a possible source of magnetic interference – was also created, but could not be added to the model for the complete aircraft. The magnetic interference was mathematically compensated for, removing the static interference, but not the in-terference during manoeuvres. The jet engine part of the complete aircraft model did not seem to generate a significant amount of magnetic interference compared to the airframe. An electric dipole, representing a submarine, was then added to the simulation. The data from that simulation was put through the mathemat-ical model and distortions of a few nT were noticeable during straight courses. The jet engine model that included rotating machinery yielded different results compared to the jet engine model in the complete aircraft model. They seemed to contain signals of higher frequency, which were however not detected by a fre-quency domain study or present during straight courses. It was thus concluded that using this particular engine model the submarine could probably still be de-tected if the course of the aircraft was kept straight, though further research is needed with more advanced models for the engine, in particular with regards to the rotating machinery.

**Acknowledgments**

This thesis is the conclusion of my studies at Linköping University. It has been a fantastic journey and I am indebted to every person I have met during my time there that has in some way made my life more entertaining and helped me de-velop from an appropriately innocent boy from Bettna to become who I am today, for better or for worse. I hope they know who they are, there are too many of them to mention.

With regards to this thesis in particular, there are a few people without whom the completion of it would have been impossible and thus deserves some recog-nition. First of all, my supervisor Alexander Hutchings, for all the assistance he has provided and all of the discussions that we have had, serious and otherwise, I give him my thanks. I want to give thanks to Torgny Hansson as well for always showing great support and interest in my work. I also want to give my thanks to the both of them for giving me the opportunity to write this thesis, which has been frustrating and humbling at times but on the whole a tremendously enter-taining, and educational experience for me. Next, I want to thank my academic supervisor, Peter Münger, for his thoughtful critique, useful comments and inspi-rational willingness to maximise the quality of this thesis. I also want to thank Tobias Hansson, my examiner, for his help whenever Peter has not been available and, of course, for his input.

There are also people that have made contributions here and there that deserve some credit. Therefore, I want to recognise the contributions of Magnus Dolk and Mats Folkesson, who have continuously provided valuable information about the system and suggestions for how to proceed. In addition to this, I want to give my thanks to Sten Rikte, Ola Borgqvist and Ola Grankvist for giving me ideas about the geomagnetic field and the submarine model. I would like to acknowledge Lennart Hasselgren as well for getting me on track with regards to the engine models.

On a note not so much related to any content in this thesis, I want to thank Hans-Göran Puke for his good company and always providing entertaining discussions. Thanks also go out to everyone else at the office that I have met during the writing of this thesis for their company, conversations and helpful suggestions. I wish them all the best of luck in the future.

Last – but most important – of all, I want to thank my family, in Sweden and otherwise, for their constant support and for laying the groundwork that made it possible for me to study and complete this program.

*Linköping, May 2019*
*Erik Olof Sandlund*

**Contents**

Notation xiii
1 Introduction 1
1.1 Background . . . 1
1.2 Purpose . . . 2
1.3 Problem formulation . . . 2
1.4 Delimitations . . . 3
1.5 Outline . . . 3
2 Theory 5
2.1 Electromagnetic fields . . . 5
2.1.1 Potentials . . . 6
2.1.2 Magnetic fields in matter . . . 6

2.1.3 Fields from a magnetic dipole . . . 8

2.1.4 Electric fields in matter and the speed of light . . . 8

2.1.5 Electromagnetic waves . . . 9

2.1.6 Electromagnetic waves in conductors . . . 9

2.2 Mathematical formulation of rotations . . . 11

2.3 Signal processing in the frequency domain . . . 12

3 Previous work 13 3.1 The geomagnetic field . . . 13

3.2 The magnetic interference of an aircraft . . . 13

3.3 Modelling the skin of the airframe . . . 15

3.4 Modelling the magnetic interference of a submarine . . . 17

4 Method 19 4.1 General . . . 19

4.2 Hardware and software specifications . . . 19

4.3 Modelling aircraft manoeuvres . . . 20

4.4 Modelling the submarine trajectory . . . 22

4.5 Modelling the generators . . . 23

4.6 Creating the mathematical compensation model . . . 23

4.7 Simulations . . . 24

4.7.1 Airframe simulations . . . 25

4.7.2 Engine simulations . . . 28

4.7.3 Main simulations . . . 31

4.7.4 Effects of rotating machinery . . . 32

4.8 Studying the frequency domain . . . 33

4.9 Approximating the electrical conductivity of air . . . 33

4.10 Sources of error . . . 34 5 Results 35 5.1 Airframe models . . . 35 5.1.1 Model 1 . . . 35 5.1.2 Model 2 . . . 36 5.1.3 Model 3 . . . 37 5.1.4 Model 4 . . . 38 5.2 Engine models . . . 38 5.2.1 Model 1 . . . 39 5.2.2 Model 2 . . . 39 5.2.3 Model 3 . . . 40 5.3 Main simulations . . . 41 5.3.1 No submarine . . . 41 5.3.2 With submarine . . . 46

5.4 Effects of rotating machinery . . . 50

5.5 Altering the electrical conductivity of air . . . 55

6 Discussion 57 6.1 Airframe model . . . 57

6.2 Engine model . . . 58

6.3 Main simulations . . . 59

6.4 Rotating machinery . . . 61

6.5 Altering the electrical conductivity of air . . . 62

7 Conclusions 63 7.1 Concluding remarks . . . 63

7.2 Future work . . . 65

A Simulation parameters 69

**Notation**

Mathematical notation
Explanation Meaning
A, a Scalar
˜
a Complex number
¯
A Vector
˜
A Matrix
i Imaginary unit
|a| Absolute value of scalar

ˆx, ˆy, ˆz Normalised base vectors

ˆn Surface normal

| ¯A| Vector norm

F Fourier transform

· Scalar product

× Product (scalars), vector product (vectors)

¯ ∇ Gradient ¯ ∇· Divergence ¯ ∇× Curl ¯ ∇2 Laplacian ∂ ∂t Partial derivative d dt Derivative xiii

Physical entities

Explanation Meaning and unit

¯

A Magnetic vector potential [Vsm−1]

¯

B Magnetic flux density [T]

¯

D Electric displacement field [Cm−2_{]}

¯

E Electric field [Vm−1_{]}

¯

H Magnetic field [Am−1]

I Electric current [A]

J Electric current density [Am−2]

¯

K Electric surface current density [Am−1_{]}

¯

M Magnetisation [Am−1_{]}

¯

P Electric polarisation [Cm−2]

V Electric scalar potential [V]

Z Electric impedance [Ω]

c Velocity of light [ms−1_{]}

d Thickness [m]

f Frequency [s−1_{]}

¯

m Magnetic dipole moment [Am2]

˜k Complex wave number [m−1]

t Time [s]
δ Skin depth [m]
Permittivity [Fm−1]
_{0} Permittivity of vacuum [Fm−1]
r Relative permittivity [1]
µ Permeability [Hm−1_{]}
µ0 Permeability of vacuum [Hm−1]
µr Relative permeability [1]

ρ Electric charge density [Cm−3]

σ Electrical conductivity [Sm−1]

θ Normalised frequency [1]

Notation xv

Abbreviations

Explanation Meaning

AC Alternating current

ASW Anti-submarine warfare

CRM Corrosion-related magnetic signature

EM Electromagnetic

EMI Electromagnetic interference

FMS Ferromagnetic signature

HSLA High strength low alloy

MAD Magnetic anomaly detector

MPA Maritime patrol aircraft

PSD Power spectral density

RMSE Root mean square error

RPM Revolutions per minute

**1**

**Introduction**

This thesis concerns the magnetic interference caused by a jet aircraft, specifi-cally with the engines placed at the rear, on a magnetometer placed nearby. This chapter begins with a short introduction and the background to the problem. It then continues to describe the purpose, aim and fundamental research questions of this work, outline the delimitations set upon it and finally give an outline for the entire document.

**1.1**

**Background**

Since World War I, submarines have been an integral part of naval warfare. The advantage of the British Royal Navy over the German Imperial Navy at the out-break of the conflict gave the latter an incentive to invest in alternative ways of fighting the British, driving forward the development of the submarine [1]. It proved quite effective: in April of 1917, almost 400 Allied and neutral merchant vessels, amounting to 900 000 tonnes, were sunk almost entirely as a result of German submarine warfare [2]. Thus, anti-submarine warfare (ASW) came to be a necessity. A multitude of resources were put into this area and a wide array of different technologies was brought forth. These were the main contributors to the dramatic decrease of tonnage of sunken vessels towards the end of the war, in particular the advances in electro-audio technologies such as the hydrophone [2].

During World War II, ASW was further developed with aircraft and radar tech-nologies enabling the detection of submerged submarines as well as surfaced ones [3]. Another technology developed during this period was the magnetic anomaly detector (MAD), which was equipped on aircraft and used to detect

merged submarines by measuring and detecting disturbances in the magnetic field generated by the Earth. During the war, these were not widely used, due to the amount of noise [4]. After the war, however, these were further developed for military and geophysical purposes, with models and technologies being pro-duced to reduce the measured noise [4, 5]. These models and technologies have been continuously applied and developed, persisting to this day in modern MAD sensors [6], which are vector magnetometers, meaning that they can measure the components and not only the intensity of the magnetic field. They also include software to compensate for the vibration of the aircraft, noise of fixed frequencies and so on [7].

The interference caused by the airframe and the manoeuvres of the airplane have been thoroughly studied, described and modelled [5, 8, 9]. Something that has not been well studied, however, is the magnetic interference from the jet engines. There are many reasons for this, but it is mainly due to the engines on most maritime patrol aircraft (MPA) being placed on the wings, far away from the MAD sensor, which is traditionally placed at the rear of the aircraft [6, 10]. An MPA currently in development, however, has the jet engines placed at the rear of the aircraft, close to a MAD of the vector magnetometer type. Thus, it has become necessary to study the magnetic influence of both the jet aircraft and its jet engines to determine whether it is even feasible to develop an MPA with its engines at the back of the airframe.

**1.2**

**Purpose**

The purpose of this thesis is to create a simulation model for the aforementioned MPA. More specifically, it aims to create a simulation model for the airframe, the engines and the magnetometer, apply the mathematical compensation model developed by Leliak [5] to the measured values at the sensor and analyse the results to distinguish whether or not the sensor system can plausibly compensate for the magnetic interference generated by the MPA. A simulation model for a submarine will then be added to determine if the MPA can feasibly detect the submarine even with the disturbance caused by the engines.

**1.3**

**Problem formulation**

With that purpose in mind, the problem formulation of this thesis can be stated as such:

1. How can the airframe be modelled? Can any simplifications be made and are they necessary?

2. How can the engines be modelled? Can any simplifications be made and are they necessary?

1.4 Delimitations 3

3. Can such a model include rotating machinery and AC generators?

4. What interfering effect do the airframe and engines have on the geomag-netic field?

5. Can the mathematical model described by Leliak [5] compensate for the magnetic interference generated by the system?

6. Given that the sensor can compensate for fixed frequencies, if rotating ma-chinery and AC generators can be added to the engine model, can a fre-quency domain analysis of the remaining magnetic interference reveal the presence of any signals related to any fixed frequency?

7. With the answers to all of the previous questions in mind: can the magne-tometer, despite its proximity to the engines, detect a submarine?

**1.4**

**Delimitations**

Several limitations are placed on this thesis. All simulations will take place in the simulation program COMSOL Multiphysics [11], inferring that no compari-son will be made with results from another software and that the method will be based on the limitations of that software. While there have been exhaustive research on the magnetic interference from airframes, with many different com-pensation models being brought forth, this thesis will solely consider the one developed by Leliak, allowing the focus to be placed on the modelling of the air-craft. Beyond this, only the magnetic interference from the airframe, engines and landing gear will be considered. From the engine, only major static components, AC generators and rotating machinery will be considered. Other components will be left to future studies. Finally, all surfaces will be considered perfect surfaces with no roughness.

**1.5**

**Outline**

Chapter 2 provides the relevant theoretical background. Chapter 3 highlights important results from similar studies.

Chapter 4 describes the method used to answer the questions as stated in sec-tion 1.3.

Chapter 5 presents the results obtained from simulation runs.

Chapter 6 gives a discussion of the results, the sources of error and the method. Chapter 7 provides a few concluding remarks and then outlines the possibilities for future research in this area.

**2**

**Theory**

This chapter gives a brief outline of the theoretical background to this thesis. It starts with describing some fundamental concepts in electromagnetic field theory and wave propagation and then goes on to mention some concepts from linear algebra and signal processing.

**2.1**

**Electromagnetic fields**

The electromagnetic field can be mathematically described as a vector field. It

has two components, the electric field ¯Eand the magnetic flux density ¯Band acts

as waves, propagating through all types of media. Central to the understanding of these waves and indeed to the entirety of electrodynamics, the field describing the behaviour of these waves, are Maxwell’s equations

¯ ∇· ¯E= ρ , ¯ ∇· ¯B= 0, ¯ ∇ × ¯E= −∂ ¯B ∂t, ¯ ∇ × ¯B= µ ¯J + µ∂ ¯E ∂t, (2.1a) (2.1b) (2.1c) (2.1d)

where ρ is the charge density, 0the permittivity of vacuum, µ0the permeability

of vacuum, ¯Jthe current density and t time. These four equations describe the

origin of and relations between the electric and magnetic fields. (2.1a) demon-strates that electric fields arise due to electric monopoles (charges), (2.1b) shows

that there are no magnetic monopoles, (2.1c) states that a time-varying magnetic field induces an electric field rotating around it and (2.1d) expresses that a cur-rent or a time-varying electric field gives rise to a magnetic field rotating around it [12].

**2.1.1**

**Potentials**

As stated by (2.1b), the divergence of the magnetic flux density is zero. Due to this, the magnetic flux density can be defined as the curl of a magnetic vector

potential ¯Aas such:

¯

B ≡ ¯∇ × ¯A. (2.2)

Inserting (2.2) into (2.1c) yields that

∇ × ¯ E+∂ A ∂t = ¯0, (2.3)

a quantity whose curl is ¯0 [12, 13]. This is then a potential field [14], meaning that it can be written as the gradient of a scalar

¯ E+∂ ¯A

∂t =− ¯∇V, (2.4)

where V is the electric scalar potential [12, 13].

**2.1.2**

**Magnetic fields in matter**

As formulated by (2.1), all magnetic effects occur due to moving electric charges – there are no magnetic monopoles. When a material is placed under an external magnetic field, the time variation of the field causes tiny current loops to appear in the material. These are very small, so small that they might macroscopically be referred to as magnetic dipoles. Now, an electric dipole consists of two charges – electric monopoles, one negative and one positive – separated by a fixed dis-tance. Likewise, if it is assumed for a moment that magnetic monopoles would exist, a magnetic dipole would consist of two magnetic monopoles – one north and one south – separated by a fixed distance. In reality, without these hypothet-ical magnetic monopoles, this corresponds to electric current loops. Magnetic materials generally consist of many magnetic dipoles. Normally, however, they are randomly ordered which usually means that they cancel each other out. How-ever, when a magnetic field is applied, these dipoles will align with it, making the material magnetically polarised – or magnetised. The resulting magnetisation is

2.1 Electromagnetic fields 7

Depending on the material, the magnetisation will align with the external mag-netic field in different ways: for paramagmag-netic materials, the magnetisation is parallel to the magnetic field and for diamagnets, it is opposite to the magnetic field. For ferromagnetic materials, the magnetisation remains even after the mag-netic field is removed. One example of such a material – from which the name is derived – is iron [12].

As stated, a magnetic dipole is in reality electric current loops. As magnetisation results in a net magnetic dipole moment for a material, it then establishes bound

currents ¯Jb = ¯∇ × ¯M in the material and bound surface currents ¯Kb = ¯M ×ˆn on

the surface. If one denotes all remaining currents as free currents ¯Jf, the total

current can be expressed as

¯

J= ¯Jb+ ¯Jf. (2.5)

Using (2.5) and assuming that the fields are static with regards to time, Ampère’s law in (2.1d) can be rewritten as

1
µ_{0}( ¯∇ × ¯B)= ¯Jf + ¯Jb = ¯Jf + ¯∇ × ¯M ⇔
⇔ ¯∇ × (1
µ0
¯
B − ¯M)= ¯Jf.
(2.6)

From this, the quantity ¯Hcan be defined as

¯

H ≡ 1

µ_{0}B − ¯¯ M[Am

−1_{],} _{(2.7)}

where ¯His the magnetic field [12].

As mentioned earlier, the magnetisation is often proportional to the magnetic field. Materials that fulfill this are called linear media. This means that the magnetisation can be expressed as

¯

M = χmH,¯ (2.8)

where χmis the magnetic susceptibility. Using (2.7), this can instead be expressed

as

¯

B= µ_{0}(1+ χm) ¯H= µ0µrH¯ = µ ¯H, (2.9)

where µ is the permeability of the material and µris the relative permeability. In

vacuum, µr = 1 holds, which is why µ0is called the permeability of free space

**2.1.3**

**Fields from a magnetic dipole**

A magnetic dipole, which corresponds to a current loop, generates a magnetic dipole moment ¯ m= ∫ S Iˆn dS= I Aˆn [Am2], (2.10)

where I is the current, S is a surface, ˆn its normal and A its area. Assuming a circular loop, this becomes

¯

m= πR2Iˆn, (2.11)

where R is the loop radius. Assuming the dipole is in the origin, the magnetic vector potential generated by the dipole in a field point ¯r becomes

¯ A(¯r)= µ0 4π ¯ m ×¯r r3 , (2.12)

where r is the magnitude of ¯r. Using (2.2), the magnetic flux density in ¯r becomes [12] ¯ B(¯r)= µ0 4π 3¯r( ¯m· ¯r) r5 − ¯ m r3 ! . (2.13)

**2.1.4**

**Electric fields in matter and the speed of light**

The magnetic phenomena described in subsection 2.1.2 have their counterparts when it comes to electric fields as well. This is, however, not of as much interest to this thesis. Therefore it can simply be noted that the electric equivalent for

magnetisation is called polarisation ¯P, currents are replaced with charges and

the electric displacement ¯Dis instead defined, which leads to

¯

D= _{0}E¯+ ¯P = _{0}(1+ χe) ¯E ←→

←→ ¯D= r0E¯= ¯E [Cm−2],

(2.14)

where χeis the electric susceptibility, 0is the permittivity of free space, ris the

relative permittivity and is the permittivity of the material [12, 13].

Also worth noting in relation to this is that the speed of light in vacuum fulfills the relation [12, 13]

c=_{√}1

0µ0

2.1 Electromagnetic fields 9

**2.1.5**

**Electromagnetic waves**

When discussing the propagation of electromagnetic fields, they are often de-scribed as electromagnetic waves. This concept can be derived from Maxwell’s equations. Assuming a medium that is non-conductive (that is, where the

electri-cal conductivity σ= 0) and with the property µr = r = 1, Maxwell’s equations

in (2.1) become [12] ¯ ∇· ¯E= 0, ¯ ∇· ¯B= 0, ¯ ∇ × ¯E= −∂ ¯B ∂t, ¯ ∇ × ¯B= µ00 ∂ ¯E ∂t . (2.16a) (2.16b) (2.16c) (2.16d) If a curl is applied to (2.16c) and (2.16d), it is found that

¯
∇( ¯∇· ¯E) − ¯∇2E¯= −µ00
∂2_{E}_{¯}
∂t2,
¯
∇( ¯∇· ¯B) − ¯∇2B¯= −µ_{0}_{0}∂
2_{B}_{¯}
∂t2 .
(2.17a)
(2.17b)
Using (2.15), (2.16a) and (2.16b), (2.17) can be written as

¯
∇2E¯= µ00
∂2_{E}_{¯}
∂t2 =
1
c2
∂2_{E}_{¯}
∂t2,
¯
∇2B¯= µ_{0}_{0}∂
2_{B}_{¯}
∂t2 =
1
c2
∂2_{B}_{¯}
∂t2.
(2.18a)
(2.18b)
These both satisfy the three-dimensional wave equation

¯ ∇2f¯= 1

c2

∂2_{f}_{¯}

∂t2, (2.19)

where c is the speed of light [12].

So electromagnetic fields in non-conductive regions are waves with a propagation speed that is equal to the velocity of light in vacuum [12].

**2.1.6**

**Electromagnetic waves in conductors**

In a conductor (i.e. where σ , 0) with dielectric properties (meaning that µ and

are not necessarily µ0 and 0, respectively), there is a free charge density ρf.

¯

Jf = σ ¯E. (2.20)

Replacing the current density ¯Jin (2.1d) with the expression for this free current

density ¯Jf yields ¯ ∇· ¯E = ρf , ¯ ∇· ¯B= 0, ¯ ∇ × ¯E = −∂ ¯B ∂t, ¯ ∇ × ¯B= µσ ¯E + µ∂ ¯E ∂t . (2.21a) (2.21b) (2.21c) (2.21d) The continuity equation for free charge is

¯

∇· ¯Jf = −

∂ρf

∂t . (2.22)

Furthermore, substituting (2.20) and (2.21a) into (2.22) yields ∂ρf

∂t =−σ( ¯∇ · ¯E)= − σρf

, (2.23)

a differential equation that has the solution

ρf(t)= ρf(0)e−(σ/)t. (2.24)

(2.24) implies that any initial free charge dissipates in a characteristic time τ =

/σ. When the charge has disappeared, when ρf = 0, (2.21) becomes

¯ ∇· ¯E= 0, ¯ ∇· ¯B= 0, ¯ ∇ × ¯E= −∂ ¯B ∂t, ¯ ∇ × ¯B= µ0σ ¯E+ µ00 ∂ ¯E ∂t. (2.25a) (2.25b) (2.25c) (2.25d) These equations have plane-wave solutions on the following form:

¯
E(z, t) = ¯E_{0}ei( ˜kz−ωt),
¯
B(z, t) = ¯B_{0}ei( ˜kz−ωt),
(2.26a)
(2.26b)

2.2 Mathematical formulation of rotations 11

where the square of the complex wave number ˜k2 = µω2+ iµσω. Dividing the

complex wave number into its real and imaginary components, postulating that ˜k = k + iκ, yields that

¯

E(z, t) = ¯E_{0}e−κzei(kz−ωt),
¯

B(z, t) = ¯B0e−κzei(kz−ωt).

(2.27a) (2.27b) This means that κ, the imaginary part of the complex wave number ˜k, affects how quickly the amplitude of the wave decreases with increasing z – a process called attenuation. The components of ˜k are

k= ω
r_{ µ}
2
r
1+_{ω +}σ 11/2,
κ = ω
r_{ µ}
2
r
1+_{ω}σ −11/2.
(2.28a)
(2.28b)
If (2.28b) is put into (2.27), it can be seen that the distance to reduce the
ampli-tude by a factor of 1/e is

d ≡ 1

κ, (2.29)

which is called the skin depth [12].

**2.2**

**Mathematical formulation of rotations**

The modelling of aircraft manoeuvres necessitates the use of mathematical ro-tations. In a three-dimensional orthonormal basis, a vector ¯v can be rotated θ radians around one of its axes ˆx, ˆy or ˆz by multiplying one of the rotation matri-ces [15] ˜ Rx(θ) = 1 0 0 0 cos(θ) −sin(θ) 0 sin(θ) cos(θ) , (2.30) ˜ Ry(θ) = cos(θ) 0 sin(θ) 0 1 0 −sin(θ) 0 cos(θ) , (2.31) or ˜ Rz(θ) = cos(θ) −sin(θ) 0 sin(θ) cos(θ) 0 0 0 1 . (2.32)

with it.

**2.3**

**Signal processing in the frequency domain**

In signal processing and analysis, an important tool for studying signals in the
frequency domain is the power spectral density (PSD). This is defined as a
mea-sure of how the power of the signal is distributed over the frequency spectrum.
The PSD can be numerically approximated with a periodogram, which is given
by
ˆ
RX[θ] = 1
N| XN[θ]|
2_{=} 1
N|F {xN[n]}|
2_{,} _{(2.33)}

where N is the number of samples, θ= f / fsis the normalised frequency, xN[n]is

the n:th sample of a signal xN and F the Fourier transform [16].

A periodogram can by its nature be very noisy. To counteract this tendency, a

method of reducing the variance of the periodograms called*averaging or Welch’s*

*method is often used. When calculating the averaged periodogram of a signal, the*
signal is divided into chunks. For every chunk, which corresponds to a frequency
interval, a periodogram is then determined. Finally, the average value of every
periodogram is calculated and used as the value of the chunk’s frequency interval
in the main periodogram [16].

**3**

**Previous work**

This chapter highlights previous studies relevant to the subject of this thesis. It begins with a description of how to model the geomagnetic field and the magnetic interference of an airframe, then goes on to demonstrate how thin layers, relevant for the simulation of the skin of the airframe, are approximated in the model and finally describes how to simulate a submarine.

**3.1**

**The geomagnetic field**

Magnetic interference, the whole basis of this thesis, depends on the fact that the Earth generates a magnetic field, the geomagnetic field, that is locally near-uniform. The distortion in this caused by conductive and dielectric materials is what is referred to as magnetic interference. The geomagnetic field in the vicinity

of Sweden has a magnitude of around 40 Am−1 with a declination of circa 71

degrees [17, 18]. It can thus be approximated as ¯

Hgeo= Hxˆx+ Hzˆz= 13 ˆx − 38ˆz Am−1, (3.1)

where ˆx is in the south-north direction and ˆz is directed from the center of the Earth [19], quite close to the one used in previous studies [8].

**3.2**

**The magnetic interference of an aircraft**

There are three primary sources of magnetic interference from an airframe:

• permanent magnetism, caused by the various ferromagnetic structural parts of the aircraft,

• induced magnetic fields created in the airframe by interaction with the geo-magnetic field,

• eddy currents, currents created in conductive components of the aircraft through the interaction with the geomagnetic field in the same fashion that currents are produced in a coil [5].

Considering these effects one by one, it can be shown that these effects result in a linear equation system with 18 degrees of freedom [5] and that the disturbed

magnetic field strength generated by the airframe itself, HD, can be formulated

as HD(t)= 18 Õ i=1 kiAi(t)= 3 Õ l=1 klAl(t)+ 9 Õ m=4 kmAm(t)+ 18 Õ n=10 knAn(t)= = HD P(t)+ HDI(t)+ HDE(t), (3.2)

where all time-dependent functions denoted with an A are functions related to projections of the field on the coordinate axes and all k are compensation

coef-ficients. The first term, HD P, corresponds to the field distortion caused by the

permanent magnetisation. The second term, HDI, represents the field distortion

caused by the induced magnetisation. The third one, HDE, characterises the field

distortion caused by the eddy currents [8].

An aircraft has three rotation angles, each corresponding to an axis of rotation: the yaw, pitch and roll angles. If these angles are respectively denoted as X, Y and Z – where all of them are 0 when they align with the geomagnetic field – the terms in (3.2) can be expressed [5, 8] as

HD P(t)= k1 k2 k3
A1(t)
A2(t)
A3(t)
= k1 k2 k3
cos(X(t))
cos(Y (t))
cos(Z(t))
, (3.3)
HDI(t)= k4 k5 k6 k7 k8 k9
A_{4}(t)
A_{5}(t)
A_{6}(t)
A_{7}(t)
A8(t)
A9(t)
= k4 k5 k6 k7 k8 k9
cos2_{(X(t))}
cos(X(t)) cos(Y (t))
cos(X(t)) cos(Z(t))
cos2_{(Y (t))}
cos(Y (t)) cos(Z(t))
cos2_{(Z(t))}
,
(3.4)

3.3 Modelling the skin of the airframe 15
HDE(t)= k10 k11 k12 k13 k14 k15 k16 k17 k18
A_{10}(t)
A11(t)
A_{12}(t)
A_{13}(t)
A14(t)
A_{15}(t)
A_{16}(t)
A17(t)
A18(t)
= k10 k11 k12 k13 k14 k15 k16 k17 k18
cos(X(t))d
dtcos(X(t))
cos(X(t))d
dtcos(Y (t))
cos(X(t))d
dtcos(Z(t))
cos(Y (t))d
dtcos(X(t))
cos(Y (t))d
dtcos(Y (t))
cos(Y (t))d
dtcos(Z(t))
cos(Z(t))d
dtcos(X(t))
cos(Z(t))d
dtcos(Y (t))
cos(Z(t))dtd cos(Z(t))
.
(3.5)

**3.3**

**Modelling the skin of the airframe**

The skin of the airframe consists of layers of mainly aluminium, graphite and fibreglass [20] that are very thin, almost exclusively between 1.6 and 3 mm [21]. The question then arises as to how to model thin layers of conductors.

Consider a slab of thickness d that is part of a geometry consisting of a material with the properties , µ and σ. Assuming that

• the material is a good enough conductor for the waves to propagate perpen-dicular to the surface alone, meaning that σ >> ω,

• the thickness of the layer is much smaller than the general geometry, • the curvature of the surface is much greater than the layer thickness, the electric field can locally be assumed to propagate only in the direction of

the normal ˆn to the surface layer, that is, if for instance ˆn= ˆx, it fulfills the

one-dimensional wave equation [22]
∂2_{E}_{¯}
∂x2 = −k
2_{E, k}_{¯} _{= ω}
r
( ± σ
iω)µ. (3.6)

If the more general case is considered, where the slab is instead multilayered, the field components of the incident wave can, when the slab occupies the region of space where x ∈ [0, d], be described [23] by the equation

Ey(x, ω)|x=0 Hz(x, ω)|x=0 = ˜Φ(ω) Ey(x, ω)|x=d Hz(x, ω)|x=d , (3.7)

where ˜Φ(ω) is the transmission matrix, which can be described [23] as a product of the transmission matrices of m layers

˜ Φ(ω) = m Ö j=1 ˜ Φj(ω), (3.8)

where for each layer the transmission matrix is given by [23, 24]
˜
Φj(ω) =
_{˜}
Φ_{j,11}(ω) Φ˜_{j,12}(ω)
˜
Φ_{j,21}(ω) Φ˜_{j,22}(ω)
= cosh(γjdj) ηjsinh(γjdi)
1

ηjsinh(γjdi) cosh(γjdj)

, (3.9)

which come from solving Maxwell’s equations for long coaxial cables and

cylin-drical shields [25]. ηjand γj, which are respectively called the instrinsic impedance

and propagation constant, are given by

ηj =
v
u
t_{}
j
µj
1+_{iω}σj
j
!
,
γj = iω
v
u
t
µjj 1+
σj
iωj
!
.
(3.10a)
(3.10b)

The surface impedance Zs, the relation between the tangential electric field at the

incident surface and the surface current density [26], and the transfer impedance

Zt, the relation between the transmitted tangential electric field and the surface

current density [26], can then be expressed as [23] Zs= ˜Φ−211(ω) ˜Φ11(ω),

Zt= ˜Φ−_{21}1(ω).

(3.11a) (3.11b) Assuming only one layer the expressions in (3.11) for the surface and transfer impedance can then be written [23, 26] as

Zs = −iωσ k 1 tan(k d), Zt= −iωσ k 1 sin(k d), (3.12a) (3.12b)

where here k= ωp( + _{iω}σ)µ.

Now, assuming that the electrical conductivity is high and that the radii of the conducting surface curvature are much larger than the skin depth δ, the Leon-tovich boundary condition [27]

3.4 Modelling the magnetic interference of a submarine 17

¯

Et= Zsˆn × ¯Ht, (3.13)

can be assumed and the tangential electric field can be related to the tangential magnetic field on either side – in the equations they are denoted by the index 1 and 2, respectively – with the relation [22]

_{¯}
Et1
¯
Et2
= Zs −Zt
Zt −Zs
ˆn × ¯Ht1
ˆn × ¯Ht2
. (3.14)

The negative sign on the impedances in the matrix comes from the fact that the

normals on either side are related to each other by a factor of −1, that is, ˆn1= −ˆn2.

(3.14) yields [28] the relations

ˆn × ¯H1=
ZsE¯t1− ZtE¯t2
Z_{s}2− Z_{t}2 ,
ˆn × ¯H2=
ZsE¯t2− ZtE¯t1
Z_{s}2− Z_{t}2 ,
(3.15a)
(3.15b)

which, using the definition of a surface impedance [13]

Zs=

| ¯Et|

| ¯K | (3.16)

and (3.13) to exchange the left-hand side for a surface current density ¯K, results

in ¯ K1= ZsE¯t1− ZtE¯t2 Z2s− Zt2 , ¯ K2= ZsE¯t2− ZtE¯t1 Z2s− Zt2 , (3.17a) (3.17b)

where the impedances are given by (3.11).

**3.4**

**Modelling the magnetic interference of a**

**submarine**

There are several contributors to a submarine’s magnetic signature. These in-clude induced eddy-currents, the ferromagnetic signature FMS and the corrosion-related magnetic signature CRM [29, 30]. The FMS is due to the ferromagnetic hull being magnetised. The CRM arises due to the moving metal hull being in electrical contact with the sea water, creating a galvanic potential difference and

inducing currents in the sea water, oxidising the metal [29, 31]. Due to this, sub-marines are usually equipped with corrosion-protection systems, which create a potential difference, called the underwater electric potential (UEP), between the ends of the submarine, yielding currents usually of a few A in size [30, 32]. There are several approaches to analysing the magnetic signature of a submarine. One of the more simple ones, that mainly considers the CRM, is to approximate the submarine with a finite current line, creating an electric dipole [31, 33, 34, 35] with a dipole length of about 100 Am [30]. Only considering the CRM can be seen as a reasonable approach, due to the possibility of compensating for the FMS via degaussing and the fact that the CRM – especially the vertically aligned one – is regarded as being among the most difficult things to compensate for [30], with protection systems leading to the UEP described previously [32].

**4**

**Method**

This chapter describes the method used to answer the questions described in sec-tion 1.3. It runs through the hardware specificasec-tions, the simulasec-tion software, how the different simulations were performed and what models were used for them, how the mathematical compensation model was produced and how alter-ing the value of the electrical conductivity of air affects the simulation result.

**4.1**

**General**

The problem of the magnetic interference is extremely difficult to calculate analyt-ically, due to the complexity of the problem and multitude of variables involved in generating it. Therefore, the magnetic interference caused by the aircraft was simulated by creating a model of the aircraft in a simulation software and placing it in a simulation sphere with the geomagnetic field from section 3.1 applied on its boundaries.

**4.2**

**Hardware and software specifications**

For all of the simulations, 32 cores on the Sigma cluster at NSC in Linköping [36]
were used. The simulations themselves were run using the simulation software
COMSOL Multiphysics 5.4 [11]. These were time-dependent studies using the
*Magnetic Fields physics interface, with the exception of the simulations *

contain-ing rotatcontain-ing parts, which utilised the*Rotating Machinery, Magnetic Fields *

inter-face. These interfaces solve Maxwell’s equations, see chapter 2, subject to certain boundary conditions [37]. The simulations were run with the standard settings,

with the exception of choosing a linear discretisation of the magnetic potential and electing the MUMPS direct solver [38]. These choices were made to speed up computation time [39].

**4.3**

**Modelling aircraft manoeuvres**

The simulations were body-centric – centered around the aircraft – meaning that any manoeuvres of the aircraft in the simulation were instead represented by rotations of the Earth’s magnetic field. The axes of the simulation’s Cartesian co-ordinate system were then aligned with the axes of the aircraft. This meant that a pitch manoeuvre corresponded to rotations around the y-axis, a roll manoeuvre corresponded to rotations around the x-axis and a yaw manoeuvre corresponded to rotations around the z-axis. Therefore, mathematically each manoeuvre corre-sponded to applying a rotation matrix from section 2.2 on the geomagnetic field in (3.1). Since the geomagnetic field was assumed to be homogeneous through-out this thesis, the surrounding magnetic field in the simulation was only altered during manoeuvres and not at all when the aircraft travelled in a straight course. In the simulations, the trajectory of the aircraft corresponded roughly to a square,

dividing it into 4 parts with a travel time of T = 7 s each – yielding a total

simu-lation time of 28 s. During each of these parts, the aircraft performed sinusoidal manoeuvres around its axes, one at a time. The argument θ in the rotation matrix relevant to that manoeuvre was set to

θ(t) = π 18sin 2πt T /14 = π 18sin 28πt T = π 18sin (4πt) . (4.1)

That is, θ was set to be a sinusoid function with an amplitude of _{18}π radians – 10°

– and a period that was equal to T/14= 0.5 s, see Figure 4.1.

The period of θ was set to T/14 due to the each quarter of the simulation, or each side of the trajectory square, being divided into 7 parts, each 1 s long. These were, in order and for n ∈ {0, 1, 2, 3}:

1. 0+ 7n ≤ t < 1 + 7n: straight course. (3.1) used for the magnetic field.

2. 1+ 7n ≤ t < 2 + 7n: yaw manoeuvre, rotation around the z-axis. (2.32) was

applied on the (3.1).

3. 2+ 7n ≤ t < 3 + 7n: straight course. (3.1) used for the magnetic field.

4. 3+ 7n ≤ t < 4 + 7n: pitch manoeuvre, rotation around the y-axis. (2.31) was

applied on the (3.1).

5. 4+ 7n ≤ t < 5 + 7n: straight course. (3.1) used for the magnetic field.

6. 5+ 7n ≤ t < 6 + 7n: roll manoeuvre, rotation around the x-axis. (2.30) was

applied on the (3.1).

4.3 Modelling aircraft manoeuvres 21 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 t [s] -0.2 -0.15 -0.1 -0.05 0 0.05 0.1 0.15 0.2 [rad]

**Aircraft manoeuvre angle**

Figure 4.1:The angle that the geomagnetic field is rotated with during

air-craft manoeuvres.

For an illustration of this, see Figure 4.2. Thus, each separate aerial manoeuvre corresponded to exactly 2 periods of sinusoidal motions as seen in Figure 4.1. When it came to the aircraft changing bearing – that is, when the aircraft had completed a quarter of the previously mentioned square – that motion was sim-ulated as a harmonic manoeuvre as well. Each turn was assumed to be a turn of

π

2 radians to the left. The turns started at 6.5, 13.5 and 20.5 seconds and ended at

7.5, 14.5 and 21.5 seconds respectively. The bearing during a turn was thus given
by the function
φ(t) = π
4(sin(π(t − 7n)) + 1) +
(n −1)π
2 *,*
n ∈ {1, 2, 3}*, t ∈* 7n − 0.5, 7n+ 0.5. (4.2)

The resulting radian from this function was then plugged into (2.32), giving a rotation around the z-axis which was applied on the magnetic field. The resulting angle between the cardinal axes of the aircraft and the geomagnetic field can be seen in Figure 4.2 below.

0 5 10 15 20 25 30 t [s] -4 -2 0 2 4 6 8

Angle between axis and geomagnetic field [rad]

**Aircraft bearing**

X Y Z

Figure 4.2:The angle between the main axes of the aircraft and the

geomag-netic field during a complete simulation.

**4.4**

**Modelling the submarine trajectory**

The submarine was simulated as being stationary. Its trajectory relative to the aircraft was thus calculated from the transformations described in section 4.3. Furthermore, the aircraft was assumed to travel with a speed of 92 m/s [6] in the positive x-direction at all times (since the x-axis was aligned with the aircraft bearing, see section 4.3), inferring that the trajectory of the submarine could be calculated as an addition of these. The simulations involving the submarine were only 7 seconds long, that is, 0 ≤ t < 7. For justification, see subsection 4.7.3. The velocity of the submarine could then be expressed as

¯vsub=

v_{0}(cos(θ) ˆx + sin(θ) ˆy), 1 ≤ t < 2
v0(cos(θ) ˆx − sin(θ)ˆz), 3 ≤ t < 4

v_{0}ˆx, otherwise,

(4.3)

where v0= −92 m/s. With a time step of ∆t = 1/ fs, the position of the submarine

4.5 Modelling the generators 23
x_{k+1}
y_{k+1}
zk+1
1
=
1 0 0 ∆tvx(k∆t)
0 1 0 ∆tvy(k∆t)
0 0 1 ∆tvz(k∆t)
0 0 0 1
xk
yk
zk
1
, k ∈ {0, 1, 2, ..., 7
∆t}. (4.4)

The initial position of the submarine was set to be at the end of the simulation sphere, 5 m from its edge, dead ahead of the aircraft and 120 m below it. The submarine entered the simulation sphere after 0.31 s, passed the aircraft after about 3 s and left it after 5.58 s. The position for all values of k was calculated using MATLAB and then imported to create a piecewise function in COMSOL Multiphysics.

**4.5**

**Modelling the generators**

The engine contained 115/200 V three-phase generators with a maximum load of 4 kW [20, 40] and an operating frequency between 324 and 596 Hz placed on the engine gearbox [20]. These were modelled by approximating the electromag-netic leakage from them with magelectromag-netic point dipoles. The size of the magelectromag-netic dipole moment was set to generate a magnetic flux density that corresponded to the limits that are set upon radiated magnetic field emissions from electrical de-vices by the United States Department of Defense standard MIL-STD-461G. The limit set upon equipment on Navy aircraft was chosen, which for 450 Hz – only 10 Hz from the middle of the operating frequency interval of the generators – lay around 114 dBpT, or 501.2 nT, 7 cm away from the source [41]. Thus, the norm of (2.13) was calculated, replacing ¯r with 0.07 ˆn, due to the amplitude of the magnetic flux density being maximised in that direction. (2.13) then became

B= µ0
4π
2| ¯m|
0.073 =
µ_{0}
4π
2m
0.073, (4.5)

which could be rewritten as

m= 2π0.07
3_{B}
µ_{0} =
343
2 ×501.2 × 10
−6_{≈}_{8.596 × 10}−2_{Am}2_{.} _{(4.6)}

Therefore, the magnetic dipole moment of the point dipoles representing the gen-erators were set to have a dipole moment of the magnitude given by (4.6).

**4.6**

**Creating the mathematical compensation model**

The mathematical compensation model for the aircraft was created by loading the simulation results into MATLAB together with vectors containing the pitch,

yaw and roll angles for the aircraft relative the geomagnetic field. These were then put into a linear least squares solver in MATLAB to minimise the function

V=Õ t | ¯H(t)| − 18 Õ i=1 kiAi(t) 2 (4.7)

for the Ai defined in (3.3)–(3.5). The compensation model could then be created

from the values of ki that gave the smallest value for V. Since these parameters

are all related to the different sources of magnetic interference, see section 3.2, their contribution could then be determined.

**4.7**

**Simulations**

The simulations performed in the software described in section 4.2 were executed in two major stages: the chronologically second of these was the one that is the main focus of this thesis. The chronologically first of these were performed to determine what models were to be put into that second stage of simulations. The focus of the chronologically first one lay on the airframe and the engine, see sub-section 4.7.1 and subsub-section 4.7.2, respectively. For description of the main sim-ulation, see subsection 4.7.3.

The properties of all materials were chosen as their standard value in the COM-SOL material library, unless stated otherwise. For all electromagnetic properties assigned to the simulation domains, see Appendix A.

The*Magnetic Fields interface in COMSOL required the electrical conductivity σ*

of all domains to be non-zero [37]. Furthermore, if the difference in electrical conductivity between two regions was too great, the computation time increased drastically. By default, COMSOL sets the electrical conductivity of air to 0 S/m

[37]. In reality its value lies around 10−14S/m [42, 43, 44]. To speed up

computa-tion time in the software, the electrical conductivity of air was set to 10−7S/m.

Water was not considered in the simulation, due to the fact that moving the body of water relative to the aircraft during manoeuvres would require implementing moving meshes, which would increase the model complexity and computation time dramatically [45]. Therefore, all surroundings were simulated as air. Thin layers, such as the skin of the airframe, were represented as boundaries

con-sisting of their respective materials. The*Transition Boundary Condition was then*

added to these boundaries, which utilises (3.17) and (3.12) to model geometri-cally but not electrigeometri-cally thin surfaces [37].

In the *Magnetic Fields interface, an initial value must be set for the magnetic*

vector potential. This was set to fulfill the relation described by (2.2). The value that was chosen was

¯
A0= _{µ}1

0

( ¯Hgeo· ˆz)x ˆy+ ( ¯Hgeo· ˆx)y ˆz =

1

4.7 Simulations 25

Figure 4.3:The first model for the airframe, as seen in COMSOL.

which fulfills (2.2) wherever µr = 1.

**4.7.1**

**Airframe simulations**

To determine how the airframe was to be modelled in the main simulation, four different models were tested, designed to increase in complexity and accuracy. Since the airframe was aligned with the x-axis, the airframe was allowed to travel in the north-south direction – that is, the x-component of the geomagnetic field in (3.1) was reversed – to allow for greater interaction between the airframe and the geomagnetic field. This, which corresponded to one fourth of the square de-scribed in section 4.3, generated a total simulation duration of 7 seconds, with the sampling frequency set to 2 kHz, at the high end of possible sampling frequencies for the magnetometer [6]. The magnetic flux density measured at the position of the magnetometer and the simulation time from each model was then compared in order to determine which model to use for the main simulation. The skin of the airframe was always assumed to be 3 mm thick, based on the information described in section 3.3.

**Model 1**

The first model for the airframe consisted solely of aluminium and simple geo-metrical objects of roughly the same size as the real airframe, see Figure 4.3.

**Model 2**

The second model for the airframe, see Figure 4.4, once again consisted solely of aluminium, but it was made more similar to the real airframe, using measure-ments from available drawings to add more detail and realism to the fuselage, the wings and the vertical stabiliser. The horizontal stabiliser was ignored due to it mainly being constructed of graphite [20].

Figure 4.4:The second model for the airframe, as seen in COMSOL.

**Model 3**

The third model for the airframe, see Figure 4.5, was identical to the second model but added the horizontal stabilisers and consisted of the main materials that the airframe consists of – namely aluminium, graphite and fibreglass. The fibreglass was simulated using the same electromagnetic properties as air, due to it being insulating. Which material the different parts of the airframe consisted of can be seen in Figure 4.5(b)–Figure 4.5(d).

4.7 Simulations 27

((a))Meshed. ((b))Aluminium components.

((c))Graphite components. ((d))Fibreglass components.

Figure 4.5:The third model for the airframe, as seen in COMSOL.

**Model 4**

The airframe is connected to several objects, some of them ferromagnetic, that could potentially have an effect on the magnetic interference of the airframe. One such ferromagnetic class of object are the landing gears. To investigate their effect, the fourth model for the airframe added steel cylinders representing the landing gears, see Figure 4.6, to the third model. In reality, the landing gears are made from HSLA (high-strength low-alloy) steel 4130 [6]. An alloy steel is, as its name suggests, an alloy of two materials, for instance of iron and carbon. A low-alloy steel would then contain a relatively low amount of carbon. The HSLA steel 4130 was approximated by choosing the low carbon steel 1020 material in COMSOL and setting the relative permeability to 1500 [46].

Figure 4.6:The steel components for the fourth model of the airframe.

**4.7.2**

**Engine simulations**

Like the airframe, different models were tested for the jet engine that were de-signed to increase in complexity and accuracy with the geomagnetic field, total simulation duration and sampling frequency remaining the same. Likewise, the simulation time was noted and the magnetic flux density at the magnetometer was measured. The rotating machinery present in the real engine was simulated as cylindrical and stationary. The former was to ensure geometrical simplicity.

The latter was due to the fact that to simulate the rotations, the*Rotating *

*Machin-ery, Magnetic Fields interface had to be used. This could however not be combined*

with the*Transition Boundary Condition, which was used to simulate the airframe.*

Therefore, actual rotating machinery was studied by itself, see subsection 4.7.4.

**Model 1**

The first model for the jet engine, see Figure 4.7 consisted of an aluminium object representing the gearbox – see Figure 4.7(a), a graphite cylinder representing the engine cowl – see Figure 4.7(b), a titanium cylinder and cone representing the

4.7 Simulations 29

engine itself and finally a titanium cylinder and two titanium cones to represent the fan – see Figure 4.7(c). Beyond this, a 1 cm thick engine cover was added. The engine cover is in reality made of some undetermined composite material [47]. For the simulation, a blank material with a relative permittivity and permeability of 1 and an electrical conductivity of 1 S/m was added.

((a))Aluminium components. ((b))Graphite components.

((c))Titanium components. ((d))Composite components.

Figure 4.7:The first model for the jet engine, as seen in COMSOL.

**Model 2**

The second model, see Figure 4.8, made the engine cowl in subsection 4.7.2 a bit more realistic by replacing the cylinder, see Figure 4.7(b), with a cylinder and two cones, all made of graphite, see Figure 4.8(a). The gearbox remained as before, see Figure 4.8(b), while the engine cover was extended to not only consist of com-posite material as before – see Figure 4.8(c), but also of an aluminum part with a thickness of 2 cm, see Figure 4.8(d). Furthermore, the engine itself – which in the previous section only consisted of a simple titanium block, see Figure 4.7(c) – was now filled out, adding compressors, turbines, the central shaft and the exhaust

cone to give it a more realistic appearance. The compressors and turbines were represented by cylinders made up of the materials that they mainly consisted of [20, 40]. For the high-pressure compressor – see Figure 4.8(e) – this was inconel alloy 718 [48]; for the high-pressure turbine – see Figure 4.8(f) – inconel alloy 725 [49] and for the low-pressure turbine – see Figure 4.8(g) – stainless steel. Stainless steel is weakly magnetic. To cover up for the worst case, the relative permeability of the stainless steel was set to 1.05 [50]. The exhaust cone and central shaft was set to consist of titanium, see Figure 4.8(i). The engine block was then set to be a 1 cm thick titanium layer, see Figure 4.8(h).

((a)) Graphite

compo-nents. ((b))domains. Aluminium ((c))nents.Composite

compo-((d)) Aluminium boundaries.

((e)) Inconel alloy 718 components.

((f)) Inconel alloy 725 components.

((g))Steel components. ((h)) Titanium bound-aries.

((i))Titanium domains.

4.7 Simulations 31

**Model 3**

The third model for the jet engine added two magnetic point dipoles to represent the AC generators, see section 4.5, to the ends of the gearboxes, see Figure 4.9. The direction of the magnetic dipole moment was set in the negative x-direction, towards the rear of the aircraft, for the generator on the right and in the positive x-direction, towards the front of the aircraft, for the generator on the left.

Figure 4.9:The third jet engine model. The purple dots that can be seen on

the gearbox are the point dipoles.

**4.7.3**

**Main simulations**

For the main simulations, airframe model 3 and jet engine model 2 were used. The discussion and justification for this choice can be found in section 6.1 and section 6.2.

In order to determine what, if any, interference arose due to the engines, two types of aircraft models were used, the first with both the airframe and the en-gine model – called the complete model – and the second with only the airframe. For each of these two simulation models, two rounds of simulations were per-formed. The first round of simulations, which can be thought of as a calibration round, had the aircraft travel in a square as described in section 4.3. The total simulation duration of that round of simulations was thus 28 s. The data from that round of simulations was then exported into MATLAB, where the mathemat-ical compensation model described in section 4.6 was created from the data and

applied on it.

With the mathematical model completed, the second round of simulations could be performed. This round of simulations added the submarine model. To keep the model as simple as possible and avoid the effects that a full moving 3D object would have on the simulation, such as greater complexity and longer computa-tion times [45], the submarine was represented by an electric point dipole with a dipole length of 100 Am, thus mainly considering the CRM of the submarine as explained in section 3.4.

The submarine was simulated so as to be stationary relative to the Earth, thus moving with the inverse of the aircraft’s trajectory relative to the Earth. Its mo-tion was simulated by adding points to the simulamo-tion along the part of its trajec-tory – see section 4.4 – that was situated inside the simulation sphere. An electric point dipole was then added to all of these points. The dipole length of each of these points was set to be 100 Am when the submarine was supposed to be posi-tioned at that particular point and 0 Am when it was not. In order to make sure that the dipole lengths were set even if the time step varied slightly, the condition was set so that the point dipole was given a non-zero value whenever the

subma-rine’s x-position was within 1_{2}v∆t of the point’s. The orientation of the dipole

moment was also rotated with the geomagnetic field. During straight courses, the dipole moment was set to be orthogonal to the bearing of the aircraft. This was due to the fact that if the dipole length was aligned with the bearing of the aircraft, the magnetic flux density measured at the magnetometer was virtually unchanged.

The round of simulations including the submarine was set to have a lower sam-pling frequency, 200 Hz, to reduce the number of points necessary for the simula-tion. Also, to minimise the computation time, the aircraft models only traversed the first fourth of the trajectory described in section 4.3 with no change of bear-ing initiated after 6.5 s. The data from this second round of simulations was then extracted and the mathematical compensation model obtained from the previous round of simulations was applied to it to give the results.

**4.7.4**

**Effects of rotating machinery**

The rotating machinery could not be added to the main simulation model due to

rotating machinery requiring the usage of the*Rotating Machines, Magnetic Fields*

interface, which did not support the *Transition Boundary Condition described*

in section 4.7. Simulations with rotating parts present would thus necessitate modelling these thin layers as domains, which drastically prolonged computa-tion time. Therefore, a model consisting of only the engines was considered for the study of the effects of the rotating machinery. Furthermore, the complexity of the second jet engine model’s geometry combined with the previously mentioned

inability to model thin layers using the*Transition Boundary Condition made *

cre-ating a simulation model in the interface using the second jet engine model quite complicated [45]. Therefore, jet engine model 1 was chosen as the basis for the

4.8 Studying the frequency domain 33 simulation model in this study, with two major differences. The first of these was that the composite engine cover, see Figure 4.7(d), was now instead made into a cylinder consisting of air and a composite layer 1 cm thick. This was due to the fact that it was previously simulated as a boundary, like the airframe skin. The second of these was that the fan – represented by a cylinder and a cone as described in subsection 4.7.2 – was set to be a rotating domain, rotating counter-clockwise with a rotational speed of 12000 RPM or 200 Hz.

The simulation was then run and the result from it studied in the time and fre-quency domain. In the latter case, particular attention was given to the values

around θ= 0.1, corresponding to the rotation frequency of the jet engine fan.

**4.8**

**Studying the frequency domain**

As stated in section 1.1, the software for the magnetometer is capable of com-pensating for fixed frequencies. Therefore, a study of the frequency domain was necessary to better justify any assessment of the magnetometer’s capability of filtering away noise. To study the frequency domain, the raw and averaged peri-odogram as described in section 2.3 was used on the amplitude of the magnetic flux density after mathematical compensation in order to determine whether any software in the magnetometer unrelated to the Leliak model could plausibly com-pensate for the remaining magnetic interference. The number of batches to di-vide the data into when using Welch’s method was always chosen as a factor of the number of data points.

**4.9**

**Approximating the electrical conductivity of air**

As discussed in section 4.7, the electrical conductivity of air was set to be 10−7

S/m in the simulations. To investigate what effect altering the electrical conduc-tivity had on the results, simulations were run with the third airframe model, as

described in subsection 4.7.1, with the electrical conductivity of air set to be 10−7,

10−5, 10−3S/m. The simulation times from these runs were then noted and

com-pared. The root means square error (RMSE) was also calculated for each value of σ: RM SEk = 1 (T+ dt) fs Õ t q (| ¯Bk(t)| − | ¯B7(t)|)2, k ∈ {3, 5, 7}, t ∈ {0, dt, 2dt, 3dt, ..., T }, (4.9)

where each ¯Bk(t)was the magnetic flux density measured at the magnetometer

of samples. The estimated ground truth, the value that was the basis for the

comparisons, was set to 10−7to save time.

**4.10**

**Sources of error**

There are several sources of error in this thesis apart from those already men-tioned throughout this chapter. First of all, there is the general simplicity of the models compared to their real-life counterparts. This can generally be mitigated by making the models more precise, though the magnetic interference from the aircraft should reasonably in the same order of magnitude. Furthermore, even though the magnetic dipole model for the AC generator is based on a worst case scenario, it is still a model with little basis in real observations. The size of this error is unknown. Also, the submarine was considered as an electric point dipole, considering only the CRM and neglecting other contributors to the submarine’s signature, such as the FMS. However, since the CRM is the most difficult part of the submarine’s magnetic signature to get rid of, it should be a reasonable as-sumption. In addition to this, as mentioned previously, water was not included in the simulation. Now, water is a dielectric medium [51] and the transition of the electromagnetic field between the media of air and water would probably have generated some interesting effects that are not included in this thesis. In all probability, the lack of water means that the submarine is easier to detect in this simulation than it would be if it included water. Beyond this, there is the chance for computation error in COMSOL. One notable such error is that the norm of the magnetic field is not computed by the software to be precisely constant dur-ing rotations, even though analytically it is. This particular error was however upon investigation revealed to be extremely small, as the deviation was not

visi-ble when a constant value 10−4Am−1off the amplitude of the geomagnetic field

**5**

**Results**

In this chapter, all the results from the simulations are displayed. It goes through the results from the simulations testing the different models for the airframe and engines, then goes on to present the results from the main simulation, the results from the simulations of the jet engine models that included rotating machinery and lastly the study concerning the effects of altering the electrical conductivity of air.

**5.1**

**Airframe models**

This section regards the simulations of the airframes as described in subsection 4.7.1. The results in this section affect what models are used in section 5.3. For justifi-cation and discussion with regards to that, see section 6.1.

**5.1.1**

**Model 1**

The simulation time for model 1 was 51 seconds. The results with regards to the magnetic flux density measured at the magnetometer can be seen in Figure 5.1. As can be seen in Figure 5.1(a), the magnetic flux density norm for airframe model 1 was constant when no manoeuvres were performed, with a value of around 57.58 µT. Oscillations in the value of the norm occurred during manoeu-vres around each axis – see Figure 5.1(b), but were least prominent for the ones around the x-axis for 5 ≤ t < 6 and most prominent for the ones around the y-axis for 3 ≤ t < 4, where the field strength reached a peak and low of about 72.55 µT and 43.45 µT, respectively – yielding a maximal deviation of around 14.97 µT

from the field strength during the straight course.

((a))The magnetic flux density norm. ((b)) The magnetic flux density com-ponents. Red is x, green is y and blue is z.

Figure 5.1:The measurements at the magnetometer for airframe model 1.

**5.1.2**

**Model 2**

The simulation time for model 2 was 2 hours, 2 minutes and 10 seconds. The results with regards to the magnetic flux density measured at the magnetometer can be seen in Figure 5.2.

((a))The magnetic flux density norm. ((b)) The magnetic flux density com-ponents. Red is x, green is y and z is blue.