[MUSIC] So far we solved the wave equation in the acoustic case, the scalar wave equation in 1D and 2D. The scalar wave equation is descriptive of sound propagation, and what I would like to introduce now is the elastic wave equation. Which is very similar, actually it's, to some extent, equivalent to the scalar wave equation in the 1D case. And it's descriptive of basically the motion of a string. Now, before we go into the actual maths of the wave equation, and the finite difference approximation, let's do a little bit of physics. I'd like to introduce the concept of stresses and the shear modulus. And for that, I'll resort to my instrument, the guitar, and we'll just make an example. So that's basically an elastic string here, and what we have, [MUSIC] If we pull the string, we have elastic waves. Now, what we need to describe the physics of that is actually the concept of stress. So if I pull here, the string in one direction, I exert a force of course, and you see that there is an angle generated. And actually the proportionality constant between the elongation and that angle is the so-called shear modulus. And that leads to the concept of stress, the stress-strain relation, which is given here. Actually, in that simple case, the stress is proportional to the strain, and the proportionality constant is the shear modulus, mu. And so on the left hand side, stress is called sigma. And because we are in one dimensions, in general its a tensor, but in one dimensions it's just one component. Sigma is equal to mu times the space derivative of the displacement u, in that case it would be uy, if we define the strain to be in x direction it would be a translation perpendicular in the y direction. So again the stress-strain relation would be sigma is equal to mu times the partial derivative with respect to x of u, the displacement. And that's going to be one of the ingredients that we put into the wave equation that we then solve using finite differences. Okay, let's use that concept of stress and develop the elastic wave equation in one dimension. Now, the wave equation constitutes, on the left hand side, the density multiplying the second time derivative of the displacement u, which equates the divergence of the stress. The stress, as we've seen before, is simply the shear modulus times the space derivative of the displacement, which is called the strength. And we keep it as sigma in that equation. Now, if we look at the entire equation, then we see that, first of all we know that the shear modulus might actually depend on x. If you think of the string that I just showed you, maybe that's not so realistic, because we assume that the shear modulus of the string is the same, it's homogeneous. Maybe not if it's a very old string, but let's for the moment assume that it's homogeneous. Actually then we can take it out of the partial derivative sign here. And then we are left, on the right-hand side, with the second derivative with respect to space of the displacement. Now, if we take the density on the right, we divide by the density in that equation, we are left basically on the right-hand side with mu divided by rho. Multiplying the Laplace operator of the displacement of the second space derivative. And now, something very important. There's a very famous result from elasticity that, actually, the shear velocity, let's call it C, is equal to the square root of mu divided by rho. So in other words, we can replace mu divided by rho by C square. And here we go, you will probably see that this equation now is basically mathematically entirely equivalent to the acoustic wave equation. We have a second time derivative on the left hand side, we have a square of velocities multiplying the Laplace operator on the right hand side. That basically just serves to show the equivalence of these two equations, at least in the homogeneous case. But of course, what's now becoming more interesting is what happens if the shear modulus depends on space, which is normally the case in general. For example, take the Earth. Of course, the rigidity is a function of space, and that will actually lead to a special solution scheme concerning the finite difference method that I would like to show you in the next couple of minutes. So starting from the wave equation that we've just derived, actually, this is a second order equation, because it contains second derivatives, both in space and time. And actually what we'll do now, we'll turn this into a first order partial differential equation, and that has a lot of consequences as to the way we will solve this numerically. In order to do that actually, first we replace, on the left hand side, the displacement by the velocity. The velocity is the time derivative of the displacement. So partial derivative with respect to time of u, and so we add that on the left hand side, and basically that replaces the second derivative by a first derivative. On the right hand side we'll keep the stress. Let's forget about the source terms at the moment. But we add another equation, and that's actually the one for the stress, which now the partial derivative with respect to time of the stress. So that's basically a stress rate, equals, on the right hand side, mu times the partial derivative of space of v. And v itself again, is the velocity, the time derivative of u. Now, this is now a coupled system of equations and actually now not only the velocity is unknown but also the stress is unknown. So this is actually a hyperbolic system of equations and this formulation is called the velocity stress formulation of the elastic wave equation. And it's very important, because many of the numerical algorithms, computer software that is today used, for example, in seismology and many other fields, also in fluid dynamics, are actually based on such a coupled formulation. And we will now proceed to find a way to solve that using finite differences. So what we'll do now actually, we introduce a new concept, so called staggered grids. In order to do that, let's go back to the original definition of finite differences. Usually we have, on the regular grid, functions defined, let's say x, x + dx, x- dx. And a first derivative approximation can be calculated using, for example, the difference of f(x + dx)- f(x- dx) divided by 2dx. That's called a centered scheme. Another centered scheme could be if we know a function, let's say, at points f(x + dx) over 2 and f(x- dx) over 2. Then we can take those two function values, take their derivative, and divide by dx. And then the derivative is actually defined in between those points and that's where the term staggered comes from. Now, the question is actually can we develop a scheme that makes use of this grid staggering? That means maybe I do not know the function and its derivatives at the same point, but I know the function at a certain grid and I know the derivative of this function. In between the function grid, and again, this is the so-called staggered-grid formulation. Why would we do that? Well, if we have basically the same computational effort we can live with a grid spacing of dx rather than 2dx for the derivative. We certainly have gained in accuracy. And so that's a very important aspect and we now proceed to try and apply this concept to our newly derived elastic wave equation. So let's do this step by step, and don't panic. At first sight, that looks a bit complicated, I admit, but you will have time later to look at this in more detail. So let's start with basically, the first part of the equation, where we have on the left hand side the ro multiplying the first derivative of velocity. And on the right hand side, we have the divergence of the stress, which is the space derivative of the stress. Now, let's first define a discretization of the stress, and you see that here in the graphics. So let's assume for the moment, and we've done that, we have discretization in space and time. And we will know the stresses at time level j and we actually know them at a space grid that is in between full indices. And we call that basically, for example, sigma i + one half minus sigma i- one half divided by dx will basically provide us with a first derivative at time level j. It will be multiplied by 1 over rho, that's actually defined at i, because the derivative of stress, so defined, is actually defined, of course, at grid points i directly. Now, the left-hand side, we have the time derivative of the velocity field. The velocities are defined in space at grid points i. And here, we actually take the difference between the velocities at time levels j + one half and j- one half. So we are at j + one half- vij- one half divided by T is of course the first derivitive of the velocity, which is the acceleration, and this equates now the previously already introduced space derivative of the stresses. Now, as you can see in the graph, both are defined at the same point. So we see here this discretization of the space time in two dimensional form, and both derivatives are staggered, but they are defined at the same point in space time, and that's very important. So what about the stresses? This is the next equation. We have the time derivative of the stress equating mu times the space derivative of vx. Now, again, if we start with the right-hand side, we will actually have velocities that are defined at time level, j plus one-half. And both, we equate, or we take the difference of the velocities that points, i + 1- the points defined at vi, divided by dx. So that's a first derivative. That's called a forward derivative, because, along the axis denoted by i, we actually look in positive i direction. And the derivative of course is then defined as point i plus one half, in between those grid points. And that's why we're multiplying it with the shear modulus, mu, defined at i plus one half. So that covers the right hand side, and on the left hand side, we have the time derivative of the stresses. Well, they are specified at i plus one half, like the right hand side. And in time we are actually taking the difference of the stresses at j plus 1 minus the stresses at j, so in time they are defined at time step j plus one half. And again you see on that two dimensional graph that both derivatives are defined at the same point in space time. Now, if we put all this together, again, just like in all other examples we had with the finite difference from relation, we can now basically, for both equations, first for the equation for velocities then for the stresses. Extract the future and express basically the future velocities and stresses as a function of everything else that we know, either at the present time step or at another time step in the past. And basically, even though that looks very complicated, this is actually an extremely elegant scheme. And the most important concept is the grid staggering is entirely consistent. We never need the derivatives at other points other than those we've introduced, but it's always important to realize that now the velocities, and the stresses and, for example, rho and mu are defined at different points in space time. And again, that's the concept of staggered grids. So on paper, that actually looks quite confusing, right? And it will need some time to actually digest this quite complicated formula. But actually, just to remove your fear, look at this simple piece of Python code that we will later also encounter in the Jupiter notebooks, just for you to see, in the actual code, this looks much simpler. Because there, you don't have to play with plus one half, minus one half, and so forth, because you will simply define vectors that are defined at different locations. So you can see here, it's very simple. You first start with a simple finite difference calculation of the stresses. Then you extrapolate the velocities, and then you take the finite difference calculation of the velocities, and you extrapolate the stresses. So it's actually a very, very elegant, simple algorithm that's extremely powerful. And again, this kind of algorithm is used in many, many different fields, and so that's why it is important, actually, to cover that topic.
