Incoherent Shallow-Water Seabed Logging Using Numerical-Model-Based Optimization for the Prediction of Unperturbed UEP Signatures

Abstract

Underwater electric potential (UEP) ship signatures are evaluated at a test site that provides only incoherent data from electric field sensors placed on the seabed. This paper shows how these data sets can be used to extract significant seabed parameters, which are evaluated by a numerical model combined with an optimization process in the frequency range from 1 Hz up to 100 kHz. This is in contrast to the controlled electromagnetic source technology used by geophysicists, which has to track both phase and amplitude, and, therefore, cannot simply be applied to data from per se stationary systems. Using the aforementioned frequency range within COMSOL Multiphysics (v6.0), the finite-element-based simulation software, a sensitivity analysis of seabed parameter behavior within a seabed depth range from 1 m to 20 m is evaluated for the highest sensitivity for each of the measurable field amplitudes of the locally available sensors. From this, the electrical conductivities of two seabed layers and the depth of the first layer can be extracted using a combination of a genetic algorithm and the Nelder–Mead method. In summary, phase information is not required to characterize the stratified seabed, providing a cost-effective approach to obtaining effective environmental conditions for free-water signature calculations.

1. Introduction

In the context of maritime security and vessel protection, the measurement of the underwater electrical potential (UEP) signature is becoming more important with regard to the rising threats of modern naval mine technology. This signature is the result of electrolytic corrosion between different metals in an electrolytic solution as well as the corrosion protection systems used. For example, the use of protective coatings (passive) or impressed current anodes (active) contributes to the resulting current flow in the solution. Together with the electrical conductivity $\sigma$, a current density J results in a measurable electric field E known as the UEP signature. In order to obtain the ship’s lowest possible UEP signature in the open sea, measurements of the ship are taken at a local site with sensors placed close to the seabed. However, the data received is subject to distortion due to the local environmental conditions. The desired UEP signature of the vessel extracted from the measurement data should ultimately result in a location-independent signature, called the free-water signature. De-embedding the environmental setup can be achieved by a seabed logging (SBL) method using incoherent AC-excited field information.

To resolve the distortions from the local seawater-seabed environment, it is necessary to determine the local stratification parameters of seawater and especially the seabed. To determine the parameters of the surrounding seawater profile, a probe can be used which measures the conductivity, temperature and the depth, called a CTD-probe in short. The more complicated task is to isolate the seabed parameters, which are also included in the field measurements. These measurements are recorded by field sensors that are placed on the seabed and, thus, receive a non-negligible influence from it. For high accuracy and successful decomposition of the measured object from environmental conditions, the extraction of the seabed parameters at the measurement site is of significant importance. Only with knowledge of the local seabed parameters can the ship’s signature be identified precisely to estimate the risk of detection at different locations at sea.

In order to obtain information about the free-water UEP signature, different methods can be used to solve the given inverse problem. In the field of underwater communications, for example, methods for locating objects near the seabed have been researched for many years [1,2]. In recent years, an attempt to compute shallow water signatures from stationary signatures was published by a group in 2022, which is based on an analytical computational model [3]. For numerical models, a PSO optimization method was presented in 2020 to perform inversions of UEP signatures associated with the representation of vessels emitting only DC signatures [4].

As a solution to unravel the local structure of the site, seabed logging techniques are used here. Different exploration methods for these techniques have been applied by geophysicists around the world for several decades now. The controlled source electromagnetics (CSEM) method is commonly used in underwater EM seabed exploration, and utilizes a source such as a dipole [5,6,7] or a coil system [8,9] under a defined applied current which acts as a transmitter and/or receiver that transmits signals into the stratified seabed and surrounding seawater, which are subsequently detected by coupled receivers. The commonly used systems are coherent systems using the phase information in time domain (TD) analysis, in which the transit times of the emitted signals in the individual layers are evaluated [10]. As early as 1968, Chessmann, Edwards and Chave presented an overview of seabed logging by using an EM source with a transient evaluation technique and outlined the relationship between the transient field distribution and the electrical conductivities in the layers [11]. Thereby, the signals used could be stationary or frequency-dependent in a very low range from 0.1 Hz up to 10 Hz, all including phase detection. Most of the systems use TD analysis, which reaches its limits in shallow water regions. There, reflections from the water surface are superimposed upon the response from hydrocarbon (HC) deposits, making it difficult to determine the exact location. To determine the parameters at shallow seabed depths, a TD method has been applied by a Swedish group to proprietary exploration systems with a coupled transmitter and receiver, whose electric dipole transmitter-receiver system can create a model of the stratified seabed at depths of a few tens of meters [12]. For locations where such phase-detecting systems are not available or cannot be built due to limited budget and experience, seabed parameters have previously been estimated to provide an approximately reliable UEP signature measured on site. The sensors on the seabed are capable of recording the values of the electric and magnetic field components. The combination of the frequency domain method and the information obtained from the local measurement system can be used to obtain a simple incoherent seabed exploration system.

In this study, a method is presented that uses the received field amplitudes measured by locally positioned electric field sensors that are excited by a sinusoidal signal for different frequencies. The parameters of a stratified seabed at the evaluated site are determined from the received signals. Due to uncertainties in the measurement, resulting from an imprecise sensor overrun of the ship, the corresponding data are pre-processed before they are used for further investigations. Based on this pre-processed data set, an environmental analysis is performed using a frequency sweep in the defined range from 1 Hz to 100 kHz, which by default includes the frequencies emitted by a ship. This procedure is evaluated using COMSOL Multiphysics software [13] based on the finite element method (FEM). The procedure represents the forward solution of the inverse problem, and, in combination with an optimization procedure, forms a search heuristic, which is to find the best fitness value between an input signal and the simulation model. The optimizer used is a combination of a genetic algorithm (GA) and the Nelder–Mead process, which uses only field amplitude data for evaluation [14,15]. Finally, different simulated input data are used to classify the optimizer, which reaches a performance of ≤1% deviation. The method is a cost-efficient way of using the hardware at a measurement site to receive information on significant local seabed parameters.

To summarize, the main research development regarding the point of innovation is given in the following short list:

• Reducing the wavelength requires a full-wave simulation model within the specific frequency range of the research 1 Hz to 100 kHz;
• Providing an identical UEP signature overrun;
• Using incoherent field data from a given site for SBL;
• Using the simulation model and optimization algorithms for shallow water SBL;

The remainder of the paper is organized as follows: Section 2 summarizes how to measure a UEP signature using E-field sensors as hardware positioned close to the seabed. Difficulties include inaccuracies due to swell and the ship’s shifting or tilting course over the sensor array, resulting in incomparable data sets. Therefore, as a precaution for further post-processing, a method is presented that normalizes the measurements and ultimately provides corrected data sets from a fixed scenario. Section 3 follows by describing the reconstruction of the measurement site within COMSOL Multiphysics. The required modules are presented, as well as possible simplifications due to the use of a dipole excitation source. Further on, in Section 4, the main analysis of the seabed parameters under investigation is presented. The derivation of each relevant seabed parameter is presented with respect to the given field data from the electric field sensors, resulting in a single frequency value for the further optimization algorithm. Section 5 outlines the parameters used by the optimization tool. The fitness function and the final optimization process of the two optimizers used are presented with the defined fitness function and the optimization results. Finally, Section 6 presents the main findings and discusses the results of the previous sections. Additionally, an overview of future applications based on this research topic is given to summarize the topic.

2. UEP Signature Recording of a Naval Vessel

This section describes the scenario of recording the electrical fields emitted by a naval vessel, which were recorded by moving past a sensor field of seven sensors located at the Aschau range, Germany. Uncertainties due to the sea state and deviations in the ship’s course are included in the procedure and are eliminated by an adaption process to provide comparable data sets for further processing. Finally, recorded sensor overruns can be made available for further post-processing without the influence of actual weather conditions and course uncertainties.

2.1. UEP Signature Measurement

The recording of a UEP signature is arranged as follows: The object to be measured (e.g., a ship) travels at a constant speed on a predetermined course to cross the stationary field sensors orthogonally. This is known as a sensor overrun. This array contains seven sensors (in the case of Aschau) placed at a fixed distance to each other, ideally on a line close to the seabed to measure the electric field. These are referred to throughout this paper as E-field sensors.

Figure 1 shows the initial situation of a measurement site, where several E-field sensors have been placed close to the seabed. The vessel’s signature is modeled by a dipole field, which is visualized as a vector field in the horizontal plane (x-y-plane) at the vessel’s depth and in the z-plane for the underkeel position.

Starting from a large distance away from the E-field sensors where the signature is below the noise level, the ship begins to cross the array. Typical distances of a sensor overrun from a normally sized vessel are 200 m, starting from a point far away and ending at a point behind the sensors, where the signature decreases and disappears again in the noise of the measurement signal. By knowing the ship’s speed, the recorded time domain (TD) fields can be assigned to a specific position.

The signature of such a naval vessel, which is generally represented by a higher order multipole, can be simplified in the far field as a dipole. Under an applied current, the dipole represents the excitation source, used below to represent the vessel’s signature, and will be referred to as the current dipole.

The seabed parameters are, of course, best detectable by placing the source as close as possible to the seabed. In order to avoid any collision with the measurement hardware placed so close to the seabed, a safe distance to the given seawater depth is chosen, resulting in a fixed source depth of 3 m below the seabed.

2.2. Adaption Process for Measurement Data to Eliminate Uncertainties of Individual Measurement Runs

During the process of performing the measurement, the recorded UEP signature is subject to uncertainties due to the swell and the displaced or tilted course of the vessel during the overrun of the E-field sensors. In order to receive comparable measured input data that can be used in the following simulation model, a correction technique is used to compensate for poor sensor overruns and provide identical data sets from multiple overruns for further processing. To obtain comparable data sets, an adaption procedure for measured electrical fields has been developed. For the Aschau site, this method is used to provide field data for a perfect measurement scenario, where the source traverses the sensors in an orthogonal direction without any deviation.

Using MATLAB software (v9.13.0 (R2022b)), the measured data from a current dipole (DC- or AC-excited) are output and displayed as the magnitude of the electric field in the x-y-plane at the depth of the sensor in a surface plot (see Figure 2). The visualized field represents the measured field of the dipole discretized by the seven E-field sensors. However, due to the large distance between each sensor, the recorded field requires better resolution to evaluate the real source position and the misalignment of the sensors. The data set shown in Figure 2 was obtained from a towed current dipole with a distance to the water surface of 1.5 m. The object was towed at a fixed course as orthogonally as possible over the sensor array at Aschau, excited by a direct current of 3 A and an alternating current signal of 75 Hz and also 3 A.

First, the measurement data is cubically interpolated in order to be able to determine the exact overrun position using better resolution. The maximum can then be evaluated to shift the electric field to manually define a specific path for the vessel, orthogonal to the sensors. Deviations in the trajectory of the overrun will result in a rotation of the field, which is calculated by analyzing the field shift using the contour lines from the MATLAB function, taking the recorded fields from each of the sensors as its input.

In order to de-embed the real sensor position on the seabed, the time difference of the maximum of the recorded field value can be used as presented in [16]. The result is a rotation angle that can be used in a coordinate transformation for field correction. An example of a pre-processed measured overrun as described above is outlined in Figure 3.

3. Simulation of the Measuring Site and Frequency-Dependent UEP Signatures

To evaluate the seabed parameters from a measured dipole overrun, we use the COMSOL Multiphysics software environment in combination with Livelink for MATLAB to investigate the signal propagation behavior within the multilayered structure [13,17]. The simulations are performed using COMSOL Version 6.0 with the electric currents (ec) physics interface within the AC/DC module and the electromagnetic waves (emw) physics interface found in the RF module. All simulations are performed using a PC with an i7-6700K processor and 64GB of RAM. When simulating frequency fields in the context of UEP signature measurements and their environmental setup, full wave simulation physics is already required due to wave attenuation caused mainly by the electrical conductivity of the given ambient media.

3.1. Representation of a Local Environmental Setup with COMSOL Multiphysics

To mimic the measurement procedure of the UEP signature in COMSOL Multiphysics, the inhomogeneous seabed has to contain a minimum of two different seabed layers ([18], p. 229 ff.). The main parameters of the seabed which influence the E-field in the given frequency range from 1 Hz to 100 kHz, namely the electrical conductivities and seabed layer depth, are to be predicted by the model. The relative permeability, which does not have an effect on the recorded E-field values, is set to 1 in both seabed layers. Also, the relative permittivity, which is less frequency dependent within the investigated range, is fixed to 30 in the first seabed layer (due to water saturation) and 1 in the second. Together with the seawater and the water-air interface, the simulation model represents four domains with different electrical parameters in corresponding layers. All domains are surrounded by an artificial domain of 10 m thickness emulating an approximately infinite domain (in the ec physics interface) or a perfectly matched layer (PML) (in the emw physics interface).

The seawater is represented by a domain of 23 m thickness (adapted to the local conditions of the measurement site) and is followed by two seabed layers: the first is set to a thickness of 10 m as a starting value for the subsequent optimization process and the second has a vertical dimension of 50 m enclosed by the PML domain. To represent a sensor overrun of a moving vessel within the simulation model, a static current dipole is used within the software, combined with a sensor line containing multiple discrete points to record each position of the virtual overrun. Each sensor is represented by a single sensor line, which together form a virtual sensor matrix to record the electric field at each point. For further simplification, and to reduce simulation time, the symmetry of a dipole field is used in combination with the following boundary conditions to make the model numerically efficient: the air domain is replaced by an impedance boundary and, furthermore, the area is split using a PEC due to the normal components of the electric field directly at the symmetry plane. The outer boundary regions are placed as close as possible to the sensor-dipole region to reduce the simulation area. This also results in a truncated sensor line at a distance of $x=-80\mathrm{m}$. The described simulation model is shown in Figure 4 and is used for all further investigations.

3.2. UEP Signature Simulation in the Frequency Domain

In order to determine the electric field values excited by a source under a defined applied current, a simulation module has to be chosen in which the effects of the field behavior within the frequency range will be included. Due to the electrical conductivities and the dimensions of the layers in Figure 4, the field distribution produces significant losses at a much lower frequency than one would expect, so that the model of the test site has to be analyzed using two physics interfaces initially [19]. Therefore, the COMSOL model used to represent the given measurement scenario is simulated using both the ec and emw interfaces. The model presented in the previous Section 3.1 is now simulated using both of these interfaces within a frequency range from 1 Hz up to 100 kHz.

For the frequency-dependent evaluation of the seabed parameters on the electric field amplitudes, the electric field values in the given frequency range are displayed in Figure 5. For the visualization of the physical comparison from COMSOL, the sensor point that is located directly below the excitation is selected in order to present the electric field behavior over the frequency range because it has the smallest distance to the dipole source.

This sensor has only an $E_{\mathrm{x}}$-component in the case of an ideal sensor overrun, and is used for the visualization of the comparison of both simulation interfaces. To show the electric field coupling between the excitation and the environmental conditions recorded at the E-field sensors, the electric field is plotted with amplitude and phase from each COMSOL physics interface. Figure 5 displays the E-field values from both simulation interfaces, representing the field data received by the sensor line. The phase information indicates the field propagation within the simulation area.

As an indicator of the need for the full-wave equation model (emw) in the underwater region, a difference limit of the field amplitudes $\delta$ from the COMSOL physics interfaces is defined in Equation (1) and set to a 5% deviation.

$$\begin{array}{c}\hfill \delta =\left|\frac{\left|E_{\mathrm{x},\mathrm{ec}}\right|-\left|E_{\mathrm{x},\mathrm{emw}}\right|}{\left|E_{\mathrm{x},\mathrm{ec}}\right|+\left|E_{\mathrm{x},\mathrm{emw}}\right|}\right|·100(\%)\le 5\%\end{array}$$

(1)

Up to the frequency point where the phase starts to differ from $0{}^{\circ }$, the magnetic field component has an influence and, thus, the ec physics interface provides incorrect results. In the water domain, which has a relative permittivity of $81$ and an electrical conductivity of $3\mathrm{S}/\mathrm{m}$, a deviation of $5\%$ is reached at a frequency of $15\mathrm{Hz}$. Above this frequency, the emw physics interface has to be used for reliable results. The phase at the above-mentioned $5\%$ threshold and the given frequency is 1.37°. The resulting zero point at the frequency of $1050\mathrm{Hz}$ is due to an interference effect of the given geometry, and, thus, the point is neither constant at the same location, nor suitable due to a non-negligible phase, which becomes close to −90°.

4. Analysis of Input Field Data Regarding the Seabed Parameter Sensitivity

In order to discriminate the individual seabed parameters from the field data received by the sensors at the measurement site, the best data point has to be elaborated in terms of the maximum sensitivity of the field amplitudes within the frequency range investigated, as well as the positioning during the sensor overrun. Furthermore, due to the fact that a numerical simulation is necessary, this analysis may lead to a minimization of the required data set and, thus, faster simulation results.

To determine the best sensitivity points of the given electric field data from the sensors, those points with the highest sensitivity regarding the variable seabed parameters (${\sigma }_2$, ${\sigma }_3$ and $d_2$) at each sensor location are investigated. For the measurement site in Aschau, Germany, this means that seven sensors placed on the seafloor can be analyzed at corresponding sensor positions at a distance of $10\mathrm{m}$ from each other. Each field sensor provides different field information due to its distance from the source, which is ultimately relevant for reconstructing the individual seabed parameters. In the simulation model, the maximum amount of different sensor data in one dimension is limited by the amount of sensors at the measurement site. Therefore, the dipole source is positioned over the second of the seven sensors as already displayed in Figure 4. The discretization of half of the sensor line was chosen to include 201 equidistant data points.

Due to the symmetrical field distribution of the dipole, the highest information content is found in a sensor chain located on one side of the symmetry axis of the dipole. By placing the source over the second sensor, it is, therefore, possible to analyze as many different sensor readings as are available from the site to obtain information on seabed parameters. If the last sensor within the row were chosen for source positioning, the adaption process from Section 2.2 would no longer be able to find the misalignment between the source and the sensor array. Due to the symmetric signature of the dipole, the sensors to its left and right of the dipole provide identical information. Over the frequency bandwidth from $1\mathrm{Hz}$ to $100\mathrm{kHz}$, the derivative of the electric field with respect to each extracted seabed parameter, electrical conductivities (${\sigma }_2$ and ${\sigma }_3$) and seabed depth ($d_2$), is calculated for each sensor as written in Equations (5)–(7) at a fixed development point. From this point, the parameters of the seabed are called the starting point and each parameter is given the index “sp” as an indicator. A two-layered seabed structure is used for the sensitivity evaluation in Equations (2)–(4) and the parameters used are given in

Table 1. Parameters at the point of investigation using seabed parameter variation.

Parameter	Value	Increment	Value
${\sigma }_{2,\mathrm{sp}}$	$1\mathrm{S}/\mathrm{m}$	$\delta {\sigma }_2$	$1\mathrm{S}/\mathrm{m}$
${\sigma }_{3,\mathrm{sp}}$	$0.1\mathrm{S}/\mathrm{m}$	$\delta {\sigma }_3$	$0.4\mathrm{S}/\mathrm{m}$
$d_{2,\mathrm{sp}}$	$5\mathrm{m}$	$\delta d_2$	$5\mathrm{m}$

along with their incremental changes. As an example, Equation (2) describes the new electrical conductivity ${\sigma }_2$ used in the sensitivity analysis in the model, which is given by the conductivity at the starting point ${\sigma }_{\mathrm{sp}}$ plus an arbitrarily chosen conductivity delta. This forms the electrical conductivity for the second parameter point for the subsequent evaluation of the sensitivity of the $\mathbf{E}$-fields of the individual sensors. The same applies to ${\sigma }_3$ and $d_2$. Furthermore, the input sensitivity value ${\delta }_{{\sigma }_2}{E}_{\nu }$ is calculated for each sensor $\zeta$ in the model and leads to the partial derivative of each component of the E-field, which depends significantly upon the sensor location (r) and the current frequency f, as well as the individual seabed parameters at the point of origin. The evaluation of the development point was carried out within the simulation model displayed in Figure 4, representing a simplification of the measurement site at Aschau, Germany.

$$\begin{array}{cc}\hfill {\sigma }_2& ={\sigma }_{2,\mathrm{sp}}+\delta {\sigma }_2\hfill \end{array}$$

(2)

$$\begin{array}{cc}\hfill {\sigma }_3& ={\sigma }_{3,\mathrm{sp}}+\delta {\sigma }_3\hfill \end{array}$$

(3)

$$\begin{array}{cc}\hfill d_2& =d_{2,\mathrm{sp}}+\delta d_2\hfill \end{array}$$

(4)

$$\begin{array}{c}{\left\{{\partial }_{{\sigma }_2}{E}_{\nu }\right\}}_{\zeta =1,2,...,7}{\left\{\frac{\partial E_{\nu }\left(\mathbf{r},f,{\sigma }_2,{\sigma }_3,d_2\right)}{\partial {\sigma }_2}\right\}}_{\zeta }\end{array}$$

(5)

$$\begin{array}{c}{\left\{{\partial }_{{\sigma }_3}{E}_{\nu }\right\}}_{\zeta =1,2,...,7}{\left\{\frac{\partial E_{\nu }\left(\mathbf{r},f,{\sigma }_2,{\sigma }_3,d_2\right)}{\partial {\sigma }_3}\right\}}_{\zeta }\end{array}$$

(6)

$$\begin{array}{l}{\left\{{\partial }_{{\mathrm{d}}_2}{E}_{\nu }\right\}}_{\zeta =1,2,...,7}{\left\{\frac{\partial E_{\nu }\left(\mathbf{r},f,{\sigma }_2,{\sigma }_3,d_2\right)}{\partial d_2}\right\}}_{\zeta }\\ \mathrm{where}\nu \in \left[x,y,z\right]\mathrm{and}\zeta \in \left[1–7\right].\end{array}$$

(7)

For each sensor, the values of the individual field components along the sensor overrun direction are determined, showing the maximum sensitivity values with dipole overrun position and frequency. Figure 6 shows the sensitivity for sensor 3, which is directly next to the source. The evaluation of the maxima are independent of any given unit, explaining why the color bars in the figure are given with an arbitrary unit (a.u.). As can be seen, the $E_{\mathrm{y}}$ and $E_{\mathrm{z}}$ components of sensor 3, the maximum sensitivity points are not directly orthogonal to the sensor, meaning that the best evaluation point can be found at a position before the object crosses the E-field sensors.

In order to obtain as much independent information from the given sensors as possible, all three field components of a sensor should be used, which is the case for sensors 3–7, located next to the intersection of the sensor line and the object path.

By selecting sensors with all components present, it is possible to select those sensors with the best electric field amplitude ratio and good sensitivity to the seabed parameters of interest. The most balanced field amplitude distributions, considering all three electric field components with respect to the sensitivity of the seabed parameters, are, therefore, not on the orthogonal sensor axis. Furthermore, they are also not directly below the source, although the highest amplitude occurs there.

To elucidate the sensitivity maxima of all seven sensors in terms of their amplitudes, the maximum of each component of all the sensors is plotted in Figure 7. Here, the distribution of the sensitivity maxima from all seven simulated sensors is displayed over the frequency range and the overrun position. Figure 7 can be used to identify identical sensitivity maxima and find a frequency and position where as many sensors as possible have the highest sensitivity. This is done to reduce the search field and, thus, the parameter space, which ultimately saves computational effort. The narrowing and possible reduction is marked by a line in Figure 7. Looking at the best sensitivity points of the outboard sensors, it can be seen that the sensitivity maximum occurs at an overrun position that is further away from the array intersection and at a lower frequency. To capture additional information from the outer sensors, a fixed frequency of $1\mathrm{kHz}$ is defined for the best sensitivity of all sensors.

It can be seen that the maxima of the $E_{\mathrm{z}}$-component are responsible for the shift of the line to a lower single frequency point. If the information from the $E_{\mathrm{z}}$ component is neglected, a higher frequency point would be displayed. Furthermore, the best dipole overrun position for the evaluation of seabed parameters lies between $0\mathrm{m}$ and $31\mathrm{m}$ away from the sensor position.

The use of only one fixed frequency point leads to a fast simulation result and saves time when taking measurements. Even though not all sensor maxima are exactly placed at $1\mathrm{kHz}$, the visualization of the sensitivities of the individual field components depicts that the field amplitude values have a slowly decreasing sensitivity and, thus, there is still a relevant intensity in the area around the maximum point.

Therefore, the frequency point of $1\mathrm{kHz}$ is sensitive enough regarding most of the sensors to obtain individual information from the field amplitudes to extract the seabed parameters. The chosen frequency and position selection are valid for this specially chosen sourcedepth. One has to keep in mind that by changing the source depth, the sensitivity has to be recalculated in order to reevaluate the frequency and best position in relation to the sensor overrun.

5. Optimization for the Evaluation of the Seabed Parameters

The following section describes how the electrical field data from the sensors is processed and combined into a single value that can be used as the input of the optimization tool. The algorithm used is presented, followed by a test function for the evaluation of the optimization process. The given incoherent E-field values resulted in a successfully adapted optimization process that fitted a seabed parameter to a given input field.

5.1. Choosing a Parameter to Define an Optimization Function

To de-embed the parameters from a layer of the environmental setup, which influences the recorded electrical field data given by a sensor overrun, the input field data is compared with the data from the COMSOL model. To achieve this, a single value is necessary that contains information about the inhomogeneous seabed derived from the amplitudes of the electric field components, which can then be used by an optimization tool to determine the local seabed parameters. This evaluation variable makes use of the incoherent field information from the measured sensor path of the current dipole. The first step is to define a test and reference system for this purpose, containing the information of all electric field components of the seven sensors. Within the input data, the test system is defined as either the measured sensor overrun or a simulation model for further evaluation of the optimization. This data set is compared to a reference system represented by the simulation model from Section 4 with variable seabed parameters to find identical field data.

This approach is shown in Equation (8) below.

$$\begin{array}{c}\Delta {\mathrm{E}}_{\nu }=\frac{|{〈E_{\mathbf{TS},\nu }\left(f\right)〉}_p-{〈E_{\mathbf{RS},\nu }\left(f\right)〉}_p|}{{〈E_{\mathbf{TS},\nu }\left(f\right)〉}_p}\end{array}$$

(8)

$$\begin{array}{c}\mathrm{Fitness}=\frac{1}{{\sum }_{\nu =1}^3\left(\Delta {\mathrm{E}}_{\nu \left(f\right)}^2\right)+1e^{-50}}\end{array}$$

(9)

where $\nu$∈$[x,y,z]$ and TS is the test system, RS is the reference system, and p represents the discrete points along the sensor line.

The inversion of the value $\Delta \mathrm{E}$ results in a fraction, which, in the best case, has a denominator of zero if the test and reference systems are identical. To avoid an infinite fitness value, an infinitesimal value is added to the sum of the field values to prevent division by exactly zero. The fitness value given in Equation (9) can be used to compare a measured dipole overrun with the simulation model replicating the measurement site. Using the results from the sensitivity analysis from Section 4, the range of evaluation of the field data from the sensor overrun is limited and used for the fitness function at a single frequency of $1\mathrm{kHz}$. The distance between the discretization points used on the sensor lines is gradually reduced from $x=50\mathrm{m}$ down to $x=0\mathrm{m}$, a range which includes each sensitivity maximum of the sensors and maintains a distance to the first maximum at $x=31\mathrm{m}$ (see Figure 7).

5.2. Optimization Tools Used

The fitness function defined in Equation (9) is used by two optimization tools, in combination with the numerical model (see Figure 4). A genetic algorithm (GA) is combined with a local fast search method, the Nelder–Mead algorithm. This combination allows the parameter search to be performed in a time-efficient manner using a comparison between the input field data set and the numerical model. Starting with the global search method of the genetic algorithm, a convergence criterion is defined to determine the input of the Nelder–Mead algorithm, which solves the function within 45 iterations maximum. The criterion is defined in Equation (10) and analyzes the iteration process using the difference between the last two fitness values over a defined iteration number. If the GA is still at the identical fitness level after several iterations, a second criterion is applied, which limits the iteration number to a specific value.

$$\begin{array}{c}\hfill \mathrm{convergence}\mathrm{criteria}=\frac{\Delta \mathrm{fitness}\mathrm{values}(\mathrm{last},\mathrm{penultimate})}{\mathrm{iteration}\mathrm{number}}\end{array}$$

(10)

In the same manner as for the GA explained in [14], the optimizer was programmed using MATLAB software and set up with the input parameters from

Table 2. Settings for the MATLAB function of the genetic algorithm.

Parameter	Value	Formula/Description
’variableLimits’	${\sigma }_2$ = $0.1\mathrm{S}/\mathrm{m}$ - $3\mathrm{S}/\mathrm{m}$;	Search field of ${\sigma }_2$
’variableLimits’	${\sigma }_3$ = $1e^{-4}\mathrm{S}/\mathrm{m}$ - $1\mathrm{S}/\mathrm{m}$	Search field of ${\sigma }_3$
’variableLimits’	$d_2$ = $1\mathrm{m}$ - $20\mathrm{m}$	Search field of $d_2$
’genBitN’	6	Bit number for the resolution of parameters from the search space
popN	$50$	Start population
iterLimit	inf	Global iteration limit
simLimit	$5000$	Simulation limit
fitLimit	$2000$	Maximum of fitness
ProgLimit	$1e^{-4}$	Progress: $\frac{Fitvalue(i-1)-Fitvalue\left(i\right)}{Iteration\mathrm{S}}$
crossOverN	$2$	Value for the bit changes in the mutation process

. A wide range of possible seabed parameters are given for the GA starting points to avoid parameter detection close to the limits, which would lead to inaccuracies in the optimization process.

The basis for the Nelder–Mead algorithm used here is the MATLAB optimisation toolbox function fminsearch, which contains the algorithm procedure [20]. This function was changed to find the maximum value for a given fitness. Its input parameters are listed in

Table 3. Settings for the MATLAB function fminsearch. The function represents a Nelder–Mead method.

Parameter	Value	Formula/Description
MaxIter	45	Set the maximum iteration steps within a constant fitness value
TolFun	1e-3	Termination criterion of the tolerance of the fitness value
TolX	1e-3	Termination criterion of the fitness parameter tolerance

.

5.3. Optimization Process with Simulated Test Data

To evaluate the optimization process generated above, several input test data sets representing possible input fields are simulated using the COMSOL model (see Figure 4). The data is used to test the fitness function in the combined optimizer. The computational time of the process is tracked and the relative deviation between the test system and the reference system is calculated in Equation (11).

$$\begin{array}{cc}\hfill \mathrm{relative}\mathrm{deviation}& =\frac{|p_{\mathrm{optimizer} \mathrm{output}}-p_{\mathrm{test} \mathrm{system}}|}{p_{\mathrm{test} \mathrm{system}}}·100\%\hfill \\ \hfill \mathrm{where}& p\in \left[{\sigma }_2,{\sigma }_3,d_2\right].\hfill \end{array}$$

(11)

The test data generated with the COMSOL model from Figure 4 is inserted into the fitness function from Equation (9) and finally entered into the optimization process (visualized in Appendix A). This was performed with a range of possible seabed parameters, which can be extracted from the literature [21,22]. For each optimization process, the final parameters and the time required to determine them provide information about the quality of the method. In

Table 4. Optimization run of the combined optimization process. A genetic algorithm has been combined with a Nelder–Mead algorithm and the corresponding termination criteria for the determination of the seabed parameters. The optimization times and relative deviations are outlined.

Pos.	TargetParameter (${\mathit{\sigma }}_2$/ ${\mathit{\sigma }}_3$/ ${\mathit{d}}_2$)	Optimization Output (${\mathit{\sigma }}_2$/ ${\mathit{\sigma }}_3$/ ${\mathit{d}}_2$)	Relative Deviation (%)	Duration
1	$1.89$/ $0.46$/ $12.3$	$1.8901$/ $0.4629$/ $12.268$	$0.01$/ $0.65$/ $0.25$	2:43 (h:min)
2	$0.87$/ $0.027$/ $13.8$	$0.870$/ $0.02740$/ $13.80$	$0.01$/ $0.56$/ $0.00$	2:24 (h:min)
3	$1.15$/ $0.23$/ $7.8$	$1.1493$/ $0.2287$/ $7.838$	$0.06$/ $0.53$/ $0.50$	2:24 (h:min)
4	$1.7$/ $0.15$/ $6.7$	$1.698$/ $0.1509$/ $6.699$	$0.09$/ $0.61$/ $0.00$	2:26 (h:min)

, data showing the relative deviation of four test systems are displayed.

These optimization processes show a relative deviation of less than 1%. The convergence time varies from $2\mathrm{h}$ $24\min$ up to $2\mathrm{h}$ $43\min$. This variation depends on the random starting values of the genetic algorithm and the fitness value resulting from the given parameters. The compact COMSOL model shown in Figure 4 requires approximately $30\mathrm{s}$ for one simulation run.

6. Conclusions and Outlook

This paper has presented a low-cost method using only the incoherent data from an existing monitoring station to obtain seabed information in a simple manner. The use of measurement data from an existing measurement site with field sensors represents a low investment. The key of using incoherent data makes the coupling of transmitter and receiver not necessary and, thus, simple data extraction at the measurement site is possible by feeding the numerical model with measured data. The most important parameters, namely the electrical conductivities of the first and second layers of the seabed $\left({\sigma }_2,{\sigma }_3\right)$, as well as the depth of the first layer $\left(d_2\right)$, can be extracted from a measured UEP signature. This was performed using a known source compared with a simulation model, in which the local parameters were then evaluated in the frequency domain. Following the analysis of previously measured current dipole electric fields from sensor overruns, a correction method for the measurement data was presented, which detects and calculates measurement errors due to the current sea state and a skewed path deviating from an orthogonal overrun of the E-field sensors. The range of parameters was then investigated to find the best positioning between the sensors and the dipole source at a depth of $3\mathrm{m}$ above the seabed that would allow a sensor overrun without colliding with the E-field sensors. The stratified seabed was analyzed at a specific investigation point using the sensor data in order to evaluate a suitable frequency range and position for the highest sensitivity of the E-field sensors with respect to the seabed parameters. This resulted in a single frequency and includes the location (i.e., necessary evaluation points) along the overun path. The absolute electric field components forming the main information for the fitness function were combined into a single fitness value using the least squares method. Two optimization tools, i.e., a genetic algorithm and the Nelder–Mead optimization method were used to determine the seabed parameters of the given local monitoring station in comparison with a simulation model.

The most important conclusion is that neither the transient evaluation techniques used by geophysicists nor phase-coupled transmitter and receiver systems are necessary to obtain local information about the seabed. In the frequency range from $1\mathrm{Hz}$ up to several $\mathrm{kHz}$, it is possible to de-embed information from amplitude signals using only a single selected frequency (here $1\mathrm{kHz}$) in the most sensitive range. The test environmental conditions can be reproduced in a simulation model with a relative deviation of less than 1%, which provides a good basis for later measured field overruns. These introduce inaccuracies due to noise components in the measured signal as well as deviations from the ideal simulation model, which should be left as unaffected as possible by the optimization system. Combined with the adaptive methodology to neglect measurement problems due to poor ship course or weather conditions, the method can be used to provide the unknown local seabed parameters at each measurement site where some local sensors are available. Knowledge of commonly influenced local parameters becomes important as it can now be used to calculate the free-water signature of a naval vessel. This is necessary in order to be aware of the vessel’s signature on the open sea and to keep it as low as possible, thus increasing the safety of the vessel and the crew.

Problems, such as rapid reflection from the water surface, as described in the introduction, are completely decoupled by the use of transient evaluation methods, as the signals emitted by the source are no longer strong enough along the entire source-to-water-surface sensor path due to the relatively high frequency. Moreover, the method used here reduces distortion due to broadband noise by choosing only one frequency and, thus, leads to higher noise immunity.

To evaluate the resulting extraction process, test data sets were generated from a numerical simulation model that demonstrated high accuracy for output seabed parameters. After this evaluation, the simulation model can be used to verify the seabed parameters at the measurement station in Aschau, Germany. The current simulation model represents this location by means of the dimensions of the layer depths and the number of sensors, which allows simple implementation in this case. To use these methods for other measurement sites, the COMSOL model must be adapted to local conditions before parameters can be accurately predicted.