Prof. David Ketcheson | Analysis and Modeling of Solitary Waves in Non-dispersive Models

Okay, great. So let's welcome our second speaker for this afternoon, David Ketcheson. David will talk about solitari waves in non-dispersive models but with a periodic structure. So David, the floor is yours. Thank you. Thank you to the organizers. It's great to be back. This talk will be kind of a continuation of the talk that I gave at the end of the program at Cambridge. Don't worry, I won't assume that anyone remembers. And before I start, let me say this is joint work with Giovanni Russo, Tania in Italy, and Blair Schoenfried at the Freshwater in Hungary, done while Giovanni was there visiting at KAUST. So, okay, so I'm going to talk about non-dispersive systems, but of course dispersion is going to enter the picture because otherwise there wouldn't be vacuum under here. It will come in through the structure of the medium. And before I do that, I just want to remind you maybe a little bit about what we expect for hyperbolic first-order hyperbolic PDEs, nonlinear first-order hyperbolic PDEs in one dimension at least. Right. And then so here's an example, and this will be the main example that I'll work with in the talk. That's the water equation. So here h is the water depth, the distance between the surface and the bottom. u is the velocity. So a drastically simplified model compared to what we were just looking at. We've averaged the depth. And g is the gravitational force here. And if you start with any initial pulse, of course it will break up into 2 pulses and then they'll steepen. And here there's nothing. to counteract that steepening, right? Only this first-order nonlinearity. And so you get shock waves. And what I'm plotting below it is the energy or the entropy, right, which decays over time after you have the shock formation, right? And if you read a book about these kinds of PDEs, they'll all have some statement like this, right? You get these head waves and the structure of the initial data disappears as you have And there's some theoretical results that prove this behavior. This is maybe the most far-reaching and well-known one, all the way back from the '60s. They say if you have a system of 2 first-order nonlinear hyperbolic PDEs and they're genuinely nonlinear, this means that the characteristic speed is changing monotonically as you Move through a wave. That essentially means that the waves appearing in the Riemann solution are shocks or rarefactions. Then for the Cauchy problem, which is what I'll be mainly focusing on, we have this theorem 3, right? If the initial data are constant outside a compact interval, the solution decays. Here I don't even care how fast it decays, just that it's going to decay to some constant state. And I'm interested— I'm showing you this because we're going to see behavior that's not like this, right? And I want to talk about why. Okay, so I'm going to use this example. I'll mention several other examples later on. So we have water waves, and the bottom here, represented by this brown stuff, is periodic. It doesn't have to be piecewise constant or continuous like this. The same kind of behavior will be observed. Only some small details will be different. Then I'm going to solve again the shallow water equations, right? Now we have this additional term here. V is the height of the bottom. Okay, and before we go to that periodic flow over periodic bottom, here's what happens over flat bottom. It's just what we talked about before, right? So I'm just getting us into the setting that I'm going to show in a moment. So this is 1 meter deep water, a 1.5 centimeter initial wave, and it's propagating meters, and we see this shock formation. Okay, right, nothing very exciting, just an end wave spreading out. Okay, here's the same experiment, but now I have these sort of submerged periodic columns on the bottom, right? So in the, in the deep parts, it's 1 meter deep. The unperturbed depth, and in the shallow parts it's only 30 centimeters. We'll watch it propagate over a similar distance. And again, as you saw there, I let the wave propagate a little bit, and then we sort of put ourselves in the frame of reference moving along with it so we can see the details. We'll do the same thing here. So there's some oscillations behind the wave, right? But we can see that instead of Shockwave during— instead of wave breaking, we get these solitary waves, or at least traveling sort of breathers, right? That's what they look like. Totally different behavior than what, what Glenn and Max tells us we should see here. Here's just a close-up of the evolution. Showing more resolution in time, so you can kind of see how these things bounce along. There's this periodic oscillation because of their phase relative to the structure of the body. And I should say, these oscillations here, as far as we can tell, and we're quite confident that they are mathematically correct, right? These are not But I won't say much about it. Okay, so does this really fall under Riemann-Lax theorem? Of course, the answer is no, because the theorem is true and we don't see the system behaving that way. So we can take— this is a genuinely nonlinear system of 2 hyperbolic equations, but there's also this spatially varying coefficient here. So one way to think about it is we can rewrite this as a system of 3 first-order equations, right? Artificially adding this equation for the bottom, changing in time. And if you look at the hyperbolic structure of this problem, right now you have 2 genuinely nonlinear fields, but also a linearly degenerate field associated with B. And this one in particular also has noncompact support because it's changing over the whole domain. So this theorem doesn't apply. From a strict perspective, it's no surprise that we get some other kind of behavior. But I think still, if you've studied hyperbolic PDEs, this is very unintuitive. It's totally different from what's generally expected. Okay. So what I want to do in the next part of the talk is to use perturbation analysis to derive some constant coefficient PDE that predicts these solitary cases. And actually, I'm really glad that my talk was scheduled on the same day after Thomas's keynote, because I think there is some connection here to discrete systems. You can kind of think of especially the piecewise constant case, right, of this is almost a discrete system.

Right.

In that sense, what I'm doing is now trying to go to a continuum system. And of course, this can, you know, help us to understand the behavior of the model, and it can also help us to more cheaply solve it numerically. So the simulations I just showed you are very expensive. You need extremely high resolution, you know, on the order of 50 to 100 points per period of be able to get the, you know, solutions that I'm showing you, which are visually indistinguishable from the exact solution. But we can solve these constant coefficient PDEs that I'm going to find much cheaper, orders of magnitude cheaper. And the perturbation approach that I'm using here was developed in the PhD thesis of Gerardo Yara. Okay, so here's just some notation. And here, you know, just to emphasize, the bathymetry doesn't have to be piecewise constant, right? It can be some smooth function. So here I've drawn it smooth. And the main thing here— so delta will be the period of the bottom, right? And lambda will be typical wavelength of our wave. And the perturbation assumption is that delta is much smaller than lambda. So we have long waves compared to the In fact, we can use the equations that we'll get to describe waves where this ratio is only about 3 or even less than that. Oh, and H is the unperturbed depth relative to the top, capital H, whereas little h is the spatially and temporally varying depth. Okay, yeah. So here's our shallow water equations. This is not a good starting point for perturbation theory because we want to derive equations that describe the macroscopic behavior, or the behavior sort of averaged over one period of this medium. But h obviously has to be varying very rapidly over that. So we're going to get rid of h and switch to eta, the surface height. And get rid of h u, which also varies, you know, as rapidly as the bottom, and switch to the discharge. Right, and those variables we get with the system equations. So now we need to— we artificially introduce this additional spatial variable, a fast spatial variable y, which— so delta is our small parameter. And y is formally defined as x over delta, but we'll treat it as being independent from x. And so then derivative partial derivatives with respect to x become this calculus, right? This sum. And we're going to assume that we have some power series for h i q, right? Power series. And every— everything that depends on y here is going to be assumed to be periodic with period delta. And our goal is to find some evolution equations for eta and q after we average over one period in y. So they're not going to describe the very fine-scale variation solution, but the larger scale. Variation. So in particular, we expect to be able to recover these solitary waves in here. Okay, so now that we've made all these decisions, there's just a mechanical process, right? It's very tedious, and it's actually quite difficult, but in principle straightforward, right? So what we do is we substitute that power series into our equations, right? Transform the derivatives. Right, based on our fast and slow scales. And then we just collect terms at each order of delta and require that they balance each other. At the lowest order, this just tells us that the leading term in the surface should be some constant. Not very surprising there. So this tells us that the waves are going to be order delta perturbations. They're going to be small waves relative to the overall depth. By the way, one thing that I'm not assuming here is that the variation in the bottom is small, right? That's large. And you saw that it was large in the simulation.

Right.

That's right. Actually, we need it to be large to get the solitary phase relative to the overall, of the same order as the overall. Okay. All right. I'm not going to focus too much on the whole process here, right? This is something that in order to do what you can do at the lowest order, we just find that the leading order part of the discharge is zero. The leading part of the surface height is constant, right? And then the next terms are independent of y. So far we've collected terms of order delta to the minus 1, order delta to the 0. We continue this process, right? Order delta to the 1 terms. Now things start to get a little more interesting. And what we do now is we solve for— so the superscripts here are not powers, right? These are indices for y. Expansion in delta, right? What I do is I solve for the y derivative terms, which will always involve the next higher order index terms, and then I integrate over one period in y. And since everything is periodic in y, the left-hand side has to vanish, right? So we get that the average of the right-hand side over one period is 0. And then there are things like— here we have the derivative of h over h cubed integrated over one period, right? That has to vanish because h is a periodic function, the distance to the bottom surface. So there are a lot of terms that appear here that then turn out to be 0 because these average coefficients vanish. And so the next order, what we end up getting is just the wave equation, but in the averaged quantities. And this is our wave speed, unsurprisingly— or sorry, this is the wave speed squared. Right, you may know the, the wave speed for the shallow water equations is square root of gh, right? Here we get exactly that, square root of gh. This is the average of 1 over h, right? But then it appears to be up there, so it has the same size. Okay, so this is basically— we've sort of linearized and averaged at this point, and as we go to the next step, we'll get into nonlinear That makes sense. Okay, so it gets really messy now. What we do is we subtract the non-average equations from the average equations at each order, right, and then integrate those up to y, and that gives us an expression for these next index variables, right? So q2 is actually just equal to its average, so it's independent of y, but Eta 2. Now eta starts to depend on this fast variable y, and there's these funny double bracket operators that appear here. They're related to certain integrals of this periodically varying function. I'm not going to worry about trying to explain this in a short time. The main thing is this— the brackets, right, mean just the average of y, and the braces mean the Fluctuating part, so the function minus its average, and then some more complicated thing that consistently appears in these equations. Okay, so now at each order we just follow these steps, right? We pull off all the terms at the next order in delta, solve for the ones with the highest index, average over one period, integrate, and and then we're. Finally, going to combine all the equations that we've gotten at these different orders into one system of evolution equations for the average surface height and discharge. And this is just to give you an idea of the sort of technical difficulty, or you might say tediousness, of this process, right? So the number of terms that you have to deal with grows complementarily just because of how calculus of products works, right? And even To get Mathematica to do this, after a while you have to be really clever. Doing it yourself, well, at least for me, is absolutely impossible.

Yeah.

Okay, yeah, but so in principle there's a huge number of terms appearing in these equations, but most of these terms, almost all of them, turn out to be 0, or they combine with each other in ways. And so we just have on the order of like A dozen terms, right? Yeah. Okay, so we define these average quantities, right, just by averaging our power series, and then we add up the equations that we found, and what we get is these 2 equations, right? So time evolution equation for the surface that is actually really simple. This one turns out to be— it's this equation to all orders, right? I've only worked it out formally to this order. This way here. And then the other equation, the one for the momentum, turns out to be much more complicated. So here are all the terms up to, you know, order delta to the 4. And what you see here is that at order delta squared now we got a 3rd derivative in time, which is going to be dispersive, and then we get even higher order dispersive terms as we continue. So We can at least mathematically justify the appearance of dispersive effects in the equations in this way. Not, not very intuitive, I would say. Does tell us that there's dispersion there. Yeah. Okay. And there's, there's like a dozen more terms in this capital F that I'm not showing you that cost me like 2 months of my life.

Okay.

So first thing you want to do now is look at the linear stability relation for these equations, right? If we first just look at up to order delta cubed, so that includes the lowest order dispersive terms, here are the equations. And if you're an expert in this area, maybe you can already see it, but this dispersion relation blows up for small wavelengths. Right? So here's the dispersion relation, right? And for small wavenumbers, you get 2 complex conjugate roots, and so one of them is widely looked at, right? So this isn't good for numerical solutions, of course, right? It's terrible. But even if you go back to like Witter's book, he tells you what to do, or well, not exactly. He tells you that if you had some other kind of 3rd derivative, that would be better. Right? And we can use these equations, which are accurate to some order in delta, to exchange time and space derivatives. Right? Because here, you know, the leading order, we have some equivalency between time and space derivatives. So we can trade some of these time derivatives for space derivatives while keeping the same order of accuracy with respect to our small parameter. And the first thing that we, that we did was to change all the higher-order time derivatives to space derivatives. And then what we get is that it's stable for small wavenumbers, but for very large wavelengths, it blows up. This is okay. You can do numerical simulations with it as long as you don't have too much spatial resolution in them. And in a way it makes sense because you know that your— these equations that you've derived are not accurate on very small scales.

Right.

What it's telling you. But what you can do instead is what's been done here, which is to just keep one time derivative in each of these high-order dispersive terms, right, and trade the rest to space derivatives. And that way you get something that is stable for all wavelengths. So this is just showing the dispersion relations of these different, you know, TTT, XXX, or XX. So this is sort of a BBM-like version of the equations. They all approximate the original system to the same order in delta, but they have drastically different things. So the simulations that I'll show you are based on this, but of course with the full nonlinear equations. You can look for traveling wave solutions in the way that you usually do. If you just look at the order delta cubed equations, now there you have only third-order derivatives. You can integrate You can, you know, look at the phase plane. You have a homoclinic connection, and you get these solitary waves that look very, very close to a hyperbolic secant squared, right? But if you include the higher order corrections in the equations, then you get something that's closer to the actual solitary waves that we see in practice. Now, I should clarify here, because we saw these, these waves before, and they didn't have like a nice hyperbolic secant squared. They weren't even symmetric, right? They were changing their shape all the time. So what we actually look at here is we, we take one point in space and watch the wave as it goes past, right? And then we get a much smoother wave in that way. And so, let's see here, the red, the outermost wave— sorry that it's so hard to see— is the one from the simulation of the full shallow water equations. The black one is what you get with a low order, just the order delta cubed in the modulation, and the dashed blue one It's closer to the red, much closer. It's what you get with the order delta to the 5th. Seems to be converging toward the real solitary shape. So I'm going to show you some numerical simulations. And I'll just comment here that for solving the original variable-theoretic shallow water equations, we use ClawPack, which is finite volume button-off type methods with numerics and rebound solvers. And I, I've refined the grid until the, you know, changes are way too small to be seen, say, in the plot. And then we solve the homogenized equations with just a pseudo-spring. So here's just a comparison with the original problem, the sort of direct numerical simulation of the variable coefficient waves in black, and then homogenized equations of different orders in blue and yellow. And this is just showing you at the bottom for scale, right, the sort of the period of the structure that it's flowing over. Not the vertical scale, of course, isn't accurate here, right? And so you can see, right, even with these order delta to the 5th terms, right, there's clear differences between Not just on the fine scale, but, you know, overall phase error between the exact solution and these homogenized equations. But we reproduce qualitatively the behavior and to some extent, you know, quantitatively as well. You get the right number of waves, about the right shape, and you can see a big improvement, right, as you go from delta 2 to delta 1. And most of all, right, this— we can solve these homogenized equations In this case, 200 times. One interesting thing is that you can actually recover even the fine-scale variation of the solution from the homogenized equations, because you get these in the process of determining the average equations. You get equations for the variation with respect to the fast scale. And so here, if you look closely, the blue is the averaged, and then The black is the exact solution, and the red is the average with the fine-scale variation added back in. And for instance, here, right, you can see it's much closer. What does it recover? Okay, so you don't always get solitary waves in these systems, right? Just because you have variable bathymetry, obviously you couldn't, because there must be, you know, as you make that variation in bathymetry smaller and smaller, You have to get back to the system that just formed an N-wave. So what happens if you have really large amplitude waves? Well, they break, right? So here's a wave that's too big relative to the variation at the bottom, and you just get— there's kind of an N-wave here, right, with a lot of incoherence in it. And we can— we've hypothesized and confirmed with computational experiments a sort of physics-based or characteristics-based criterion for when you'll get shock formation versus when you get solitary waves. Another way to think about this is that there's kind of a limit to how big the solitary waves could be in the system.

Right.

They, like KDB solitons, they get narrower as they get taller, right? But there's an inherent spatial scale of the system, and eventually they sort of are, you know, weak reaching the limits of that spatial scale. So maybe there are some sort of solitary wave-like structures appearing out of this eventually. There's lots of other physical settings where this kind of thing can arise. So originally the paper that I read that kind of got me interested in this is this one by Levesque and Yang from 2003 where they looked at elastic waves. Really they looked at the P system, right? So But they use the notation of elasticity there. My master's student, Leila Rosala, in her master's thesis looked at the gas in a one-dimensional pipe with a spatially varying cross-section where that spatial variation is periodic. You get the same kind of behavior. We've looked at water waves where the variation is in the transverse direction, right? So water waves in a channel that's not flat bottom. Actually, what's interesting, we're going back at least to Peregrine. And then I got really interested in this paper by Sergei, who's here, and Giovanni Russo and his student, where they looked at a layered gas, right? So where you have alternating layers of 2 different gases. And they didn't see exactly the behavior that I've been showing you today, but they were looking in sort of the large, very large amplitude wave regime. If you look at smaller perturbations, then you do get this. So here's— these are the compressible Euler equations. It's another first-order system of hyperbolic PDEs, right? Generically, you get shock formation and entropy decay here. But here I'm solving— my initial pressure is constant except for this pulse right here. Here the pressure is constant, and then I have sinusoidally varying density. And also temperature, right? They compensate for each other in that constant measure. And I'm sending this pulse into this sort of oscillatory background medium, and instead of shock wave formation, we get the emergence of solitary waves because there's a backscattering distortion. This one's maybe more interesting on the surface at least because the density is also evolving in time here, right? It's not like the bottom of our shallow water equations. But we can think about this in a way that really brings it back to being very similar to the shallow water equations, and that is that if you write these equations in the Lagrangian form, and if you assume that there's no shock formation, then yeah, then you can rewrite this as the P system, right? Where here V is now the specific density, right? So your specific volume, one over the density. And S here is the entropy, which is just a function of x. And this has almost the same structure as the shallow water equations with this spatially varying, spatially periodic coefficient. And if we introduce another equation for S formally, then it's a system of 2 generally nonlinear characteristic fields and one linearly degenerate field. Really has same structure as the shallow water equations. And actually, like I mentioned, Levesque and Young looked at the P system back 20 years ago. So you can almost just directly read off the homogenization of the Euler equations from their work after some transformations and get something like this here. I've already put this in sort of the mixed space of time-relative form that we know is more stable. Again, there are some higher-order nonlinear terms here that wouldn't fit on the slide. And here's finally just a comparison So in blue and in yellow is the Euler equations of the P system, which match perfectly. They should match as long as there's no shock formation here. Confirmation of that. The dashed green line is what we get from these homogenized equations. Okay, so the main point is just that periodic variation of coefficients can induce an effective dispersion in first-order hyperbolic systems and make their solutions behave like, you know, the kind of waves that are being discussed in this room. And this can arise in a lot of different physical applications. There's a ton of open questions. One that I— yeah, I've been thinking, you know, today more, you know, are there connections here to discrete systems? Another interesting one is that there seems to be sort of a limit to how fast these waves can go. And this has to do with, right, in the first-order system, you have Right. And so there's also an amplitude limit to solitary waves, and they seem to be compact, or at least their front edge has a kink, right? So there's some interesting questions there. And of course, just generally, you know, can we prove that these are continuous solutions for all time? I'll stop there. Thank you, David. Uh, we have time for a few, few questions. Yes. Did you try to compute them as exact solutions, like periodic orbits, so to speak, something like that, periodic modular shifts, so that we could— No, I want to, and I talked to Mark Hofer about this 2 years ago, but so far I haven't succeeded. But I haven't tried particularly hard, and I don't think I know the right tools that well, but Maybe you do. I thought so. Yeah. Yeah. Thank you very much for a nice talk. So I'm aware that you use idealized equations, of course, to tackle the problem. But generally, if you have such kind of settings, you expect to be having reflections. We have instant wave effectively transmitted. So from your simulations, you know, this is actually put on by some studies. Now, have you seen any effects or any wave reflection? Oh, there's reflection going on everywhere here. Yeah, absolutely. But so you have waves traveling in both directions and the net wave that you see are these solitary waves. I haven't tried to give an intuitive explanation of that, but, but yeah, reflection— you wouldn't have this behavior without reflection at all. Yeah, you just get like burgers. I was expecting having— if you start with the wave train or the solitary wave train, you'll see really something going backwards quickly. So of course, if your wave— if its wavelength is short relative to the period of this, especially in the piecewise constant case, right? If you have a small wave compared to that, then you would see the reflections, right? Very clearly. But these waves are at least a few periods long. Right? So you don't see individual reflected waves. But, but yeah, this dispersion is really caused by reflection, 100%. Yeah, you can even— so you can, you can see it more clearly if you analyze just linear wave equations, right? Like the linear wave equation with, you know, periodic coefficients, right? It's this reflection leads to dispersion.

Thanks.

Yes.

So, you know, in your simulations on the shallow water system, You see this clear radiation behind. Yes. Which I think reminds me of, like, you know, if you have an inhomogeneous, for example, discrete system with, for example, alternating masses, or, you know, some kind of inhomogeneity, then you would actually— then there is actually no generic solutions, right? And you would see this kind of radiation coming out of them. So I wonder if that's— Yeah.

So, okay, so first of all, I didn't show any solitary wave collisions here, but they do behave roughly like soliton collisions. It's not clear whether there's radiation in that case. Here, just, you know, sort of a large hump resolving and, you know, breaking up into solitary waves. We do see a lot of oscillations. One hypothesis I have is that because especially at the beginning you see what looks almost like a standing wave, that I think it's sort of almost linear wave that's maybe in the band gap for this medium, right? So that's one hypothesis.

Maybe it's just shedding some energy.

And so maybe it's constantly shedding stuff that's in a band gap or something like that. But that's something I'd like to understand better.

The average equation, you don't see the radiation.

Yeah, those disappear in the average equations. Absolutely. Okay, final question, and then we'll go. So usually when— because you're going to high— this is the high-order equation theory. Usually if you go too high, typically this is a simple series. So at some point you just start getting worse. Yeah, I mean, I did not go too much in that photo, but I wonder if you will look at such effects. Basically, I think in the equation If that's what you put in the initial condition, then contribution from Kenneth and plus other terms becomes comparable. Then things would be wrong. You should stop. Did you see any such effects? Not yet. I don't. Yeah, so like I said, it took me like a month to work out the highest order terms that are here, so I don't anticipate ever going further than that. So yeah, I don't know if it's an emergent sequence. Probably not. But I don't know a place else to check when you should talk to each other.

Yeah.

Okay, let's— thanks. We'll speak again.