Table of Contents

## Introduction

How do you determine the shape and strength of all the major reflectors (reflectors are the transitions between rock types) deeper and deeper inside the earth? Reflect sound waves and feed the raw recordings into the world’s best algorithm for unravelling them … This is XWI (the algorithm) created by S-Cube – founded by researchers at Imperial College London.

Think of S-Cube’s avant-garde algorithm – XWI – as the * brain behind the earth’s audio. *Our human brain processes the sounds that we hear, allowing us to get a coherent signal of the direction that the sound came from, and other information (such as the words that were said). Similarly, this is analogous to the receivers in the earth, which receive a signal and pass it into XWI, which acts like the brain by processing the waves and telling us what is in the subsurface.

XWI does this in a completely unique and innovative way which means that they hold mathematical commercial advantages when pinpointing natural resources deposits (figuring out which reservoirs in the earth have oil and gas), whilst also even storing carbon (underground!) to offset the emissions. The fact that the algorithm can determine which reservoirs are suitable for storing carbon in the subsurface means that we can make the Earth much more clean for the future by reducing emissions.

During my work experience at the company S-Cube, I was set the challenge of working on a highly simplified version of a related problem to which XWI tackles (parameter extraction with 1 unknown instead of several millions in real-life!!), which I found extremely enjoyable solving in Python. My project was titled “Parameter Extraction for 1 unknown”.

In particular, I was exposed to some very important new mathematical concepts, such as the Newton-raphson method, which was thanks to my supervisors: Dr. Nikhil Shah and Dr. John Armitage, who both read Mathematics at Cambridge and Oxford respectively.

To understand the situation, in Figure 1 below we see that the source sends a sound wave down towards the seabed. Once it hits the seabed, some of the sound waves are reflected back up towards the surface where the receiver receives the waves. (*some waves pass through the seabed to the rock layers*)

In the first case (Figure 1), the distance between the source and the receiver is set to be 0, so that the depth (* h*) is simply the vertical distance down to the seabed.

In Figure 2 below, we can setup the source and the receiver at a set distance away from each other, where the * offset (x)* is the distance from the source in Figure 1 to the source in Figure 2.

## Slowness

The term ** slowness** (

*represents the*

_{S1)}**speed that the wave travels**at in the Earth.

This is a very important parameter to map because *slowness *shows us all of the geological layers, which form a trap to store a fluid. In Figure 3 below, we can see that *high slowness* is trapped between *low slowness* layers – in order to keep the fluid reserves in place (e.g. oil)

If we find out the slowness within a layer, we can tell if that layer contains more fluid (*e.g. oil*) which is useful to know.

Slowness is the **reciprocal of velocity**:

Our **aim **is to find out the **travel time** for the sound waves from the moment the source sends the wave to the moment the receiver receives the reflected wave.

Thus, we can use the following equation to help us work this out. A simple equation is the Time = distance/velocity. Since *slowness* is the reciprocal of velocity, we have an equation for travel time in terms of distance and slowness.

In this way we can work out the *travel* *time* for both Figure 1 and Figure 2, which we will complete below.

## Comparison to Image Classification

In order to understand the Parameter extraction task more easily, we will compare the task to *image classification* in Figure 4 below.

Parameter extraction (inversion) | Image Classification | |

Model | slowness, reflector (h) | Neural network |

True data | travel times | What the images are of |

Input | pulse (from the loudspeaker) | Images |

Predicted Data | model travel times | What the neural network thinks the images are of |

Updated model | change slowness | Update classification of image |

AIM | To find the values of the true model (hidden slowness value) | Classify image correctly |

## Overview of the model

The simple diagram of how the program will work is in Figure 5. We **define the true data **(which is hidden to the program, as this is what we are trying to find).

Next we **start our model **with a guess input, which outputs predicted data. We calculate the **error **(*misfit)* between our predicted data and the original true data, and then update the model accordingly.

These last 2 steps are **looped** until the predicted data is as close as it can be to the original **true data**.

## Overview of our Parameter Extraction problem

Figure 6 shows the overview of the steps taken by the program to complete the *Parameter extraction* problem. I will go into more detail for each step below.

#### Define True Model

In our example, the **true original data** is simply the observed data from a “real-world” seabed. We can work out the depth, and the slowness values – and store this into our program. (*this data is ‘hidden’ because using our model we want to get closer and closer to this value from a random guess input value for the slowness)*. Figure 7 shows us define this *true data*:

#### Calculate recording (for travel time)

The next step involves taking the true model that we defined above and working out the True travel times (* t twiddle/tilde*) for the model.

In order to work out the **travel time** for the diagram in Figure 1 (*where distance between source and receiver was 0)*, we use the equation below. The *true_depth (h)* multiplied by 2 because the total distance travelled by the sound wave is down to the seabed and back up again from the seabed.

Now that we have the distance travelled (*2d*), we can use the following equation to work out the travel time for (t twiddle x):

The code for the steps stated above is shown in Figure 9. In order to use the “math.sqrt” function, we had to import the math library at the start of our function.

#### Make our starting model

Now that the *true model* has been defined, we make our model by firstly **defining guess_slowness**, which is what our model will start with. The aim of the model is to keep looping around the model in order to get this *guess_slowness* value as close as possible to the *true_slowness* value defined above.

**Defining our model depth **(see Figure 11), and the code for it in Figure 12.

#### Calculate recording (for travel time) using our Model

Now that we use got our *guess_slowness *and *guess_depth* values, we need to calculate the travel time using these values, similar to the steps in Figure 9.

In Figure 13, “tpx” simply means the predicted travel time value for the x-offset value. “d_guess” is now our guess distance value using our *guess_depth *value which we defined above.

#### Calculate Misfit

A **misfit function** simply measures how closely our model’s predicted values are to the true values.

The formula for misfit is shown in Figure 14. Simply put, we are comparing our “tpx” value (predicted travel time) which we calculated in Figure 13 to the “t twiddle x” value (true value) which we defined in Figure 9. This will tell us how much **error** we have created compared to the true value.

#### Change guess_slowness value

After calculating this misfit, we can now change our *guess_slowness *value by a small amount (10^-2) – and then see what effect that has on the misfit. The code in Figure 15 shows how adding 0.01 to the previous guess_slowness value, and then we recalculate our travel time value.

#### Recalculate Misfit and Compare

After changing our predicted travel time value, we can now **recalculate our misfit** to see whether the change we made was beneficial or not.

Figure 16 shows the code used to compare the two *misfit* values.

If the difference between the old and new misfit value is greater than 0, it means that the misfit_new is smaller and * the smaller the misfit, the more accurate the slowness value, because it is closer to the true slowness value*. This means that the direction we want to go is positive (we can keep adding 0.01 to the guess_slowness value until we get close to the true value).

See Figure 17 to see how we must travel in the direction which will take the graph to y = 0.

If the difference between the old and new misfit is less than 0, it means that the *misfit_new* is larger than the *misfit_old, *which means we are going in the wrong direction, because we are getting further away from the true value of the slowness. This means that direction = -1, so that the program knows that we need to subtract 0.01 each time to get closer to the true value.

#### Update guess_slowness in the correct direction

Now that we have worked out the direction in which we want to go, we can simply update the direction in which we change *guess_slowness* each iteration.

Figure 18 shows the code which we use to update the direction.

**Delta **stores the **absolute value** of the difference between *misfit_old *and *misfit_new*:

To find out what we need to update, we can follow the steps shown in Figure 19:

#### Record guess_slowness in each iteration

As shown in Figure 6, there is a loop in our program. This loop will keep updating our model and checking that the misfit is getting smaller and smaller. Just for the purpose of the project, the loop I made will iterate through 10 times:

At the bottom of the loop (before iterating again), the new guess_slowness value is **appended **to an array called ** slowness_array**. This is so that we can plot the graph in the next step.

#### Plot graph of guess_slowness values

Finally, to show our model working in action, we want to **display a graph** to show our guess_slowness values getting closer and closer to the true value (where y = 0 is the ideal!).

To plot a graph we import another library called **matplotlib** at the start of our code, and then we can use the code in Figure 20 to plot the graph:

The graph that we get from our code is shown in Figure 21, which shows that over our 10 iterations, the value for guess_slowness gets more accurate each time – due to our calculations using the misfit values.

## The Newton-Raphson Method (Newton’s method)

Our project above finds the minimum of f(x) in order to get an accurate answer for the *guess_slowness *value.

The Newton-Raphson method is an algorithm used to *efficiently *and *quickly* produce successively better approximations to the **roots of a real-valued function**. The method involves using straight-line **tangents** to find the approximations of the result.

Essentially, a summary is **finding x such that f(x) = 0**.

### How does the method work? Visual representation:

The end goal is to find the point where f(x) = 0 – i.e. the point marked *g* in Figure 22.

To start with, we make an **initial guess** called x1. We calculate f(x1) and the most likely case is that f(x1) is not equal to 0 (if it is… then we are very lucky and our job is complete). This is show in Figure 23.

Now that we have tried x1, we now want to try again with a better guess. Newton’s method says that we should draw a **tangent **to the graph at the point (x1, f(x1)) for a * good linear approximation*, and where this tangent crosses the x-axis is our x2 – our next guess. See Figure 24.

Essentially, this process is continued for x3, x4, x5 … until we hit point *g*, where f(x) = 0. In figure 25, you can see that with this method, we are quickly getting closer and closer to the **root** of f(x).

### How does the method work? Mathematical representation:

The most important step in Newton’s method is to find where the tangent to x1 meets the x-axis, because this gives us our next better guess.

We need to find this x-intercept, and thus we need to find the **equation of the tangent line**. Figure 26 and 27 show how we can find the equation of the tangent using the equation of a line.

Thus, we can now find the x-intercept of the tangent, as shown in Figure 28:

### Newton’s method Equation

Therefore, by following the exact same method – we can now show the equation of Newton’s method **algebraically**, as shown in Figure 29.

## Conclusion

This project would not have been possible without the kind support, encouragement and explanations by the team at S-Cube, specifically **Dr. Nikhil Shah** and** Dr. John Armitage**. Thank you very much for your time and guidance. I learn several new skills and topics, as well as learning how to link mathematics to real-world solutions, such as S-Cube’s award-winning XWI waveform algorithm.

My project now is to solve the same parameter extraction problem, but instead for** multiple layers where slowness and depth are unknown** (instead of just slowness).This is a more realistic scenario in real-world earth surveying, because there are multiple unknowns to solve for.

For example, in the Mona Lisa image below, there are millions of pixels, where each pixel is an unknown to be found. It requires a very complex algorithm and powerful computer to solve for one million unknown simultaneously!

I hope you enjoyed reading this article – please feel free to leave any comments/feedback below (in the comments box). Thank you!