 Let's start by thanking the organisers for the opportunity to speak at this interesting meeting. It's the third one I've been to and they've always been most stimulating. So this is a collaboration between Caltech, Calfton University of Queensland. The work I'm going to describe was mostly done by my post-op voter, Mostat. So it's going to be divided into three parts. I'm going to talk about the technique known as shop dynamics. There's a method developed by G.B. Wittem in the back in the 50s to describe the dynamics of shock waves. So I'm going to give a little primer on that first for those who aren't familiar with. Then the second part going to be the application to the evolution of both planar and MHD shocks and also cylindrical symmetric shocks. And then I'm going to give the third part is going to be an analytic basis for some of the findings in part 2. So I'm going to talk about the method of shock dynamics. This is a simplified model for determining the behavior and motion of shocks in 1, 2 and 3D invented by G.B. Wittem back in the 50s and described in detail in his book in the 70s. The idea is to obtain the behavior of shock waves by considering a shock moving down a tube with the area that varies with distance. What Wittem actually did was write down the one-dimensional equations of 1D flow with area change and solve the equations linearly using the method of characteristics. Then substituted the shock dynamic conditions across a shock and obtained an equation to describe the shock muck number as the shock propagates. But he later discovered that he could do this in a more precise way. He formulated the characteristic rule which is write down the nonlinear differential equation for the C minus characteristic end during the shock from behind. Substitute the expressions from the shock jump conditions directly. The resulting differential equation gives the variation of muck number with radius. An example of this would be the schematic for an RT diagram for a converging shock. Here is the convergent part of the shock. The origin of the shock reflects back off the origin. In this case, a relevant characteristic is the one that enters the shock from behind. This is just one example of the converging shock application. Shock jump conditions for gas dynamics are given by these relationships for pressure and density. There's also one for velocity, but it's contained in the density relationship. Substitute these into the C minus equation. You get a characteristic equation. You get an ODE for m as a function of x. If a as a function of x is known, this can be integrated to give this formal solution in the strong shock limit from up number greater than one. This lambda m which appears is a known function of m. This approaches a constant which is just a function of gamma. I'll call that Wittem's constant. This ODE then can be integrated to give the well-known result for the goodly collapsing shock. Wittem's approximation gives the exponent given by this value compared with the exact value. The approximation is very good for converging shocks. For MHD, one can do the same thing, but it's more complicated. The starting point is the ODE or the PDE is for a one-dimensional flow with area variation, but the equations are more complicated because they involve magnetic fields plus equations for the dynamics of magnetic fields. You can write down the appropriate, I've got C plus. You've got to find the appropriate physical characteristic. There's a shock from behind. You can write down the left eigenvector and you can obtain the characteristic equation given by this expression. I haven't written it down in four, but it's quite straightforward, but much more complicated than for the gastronomic case. To get closure you need, according to Wittem's rule, you need the shock relations which are much more complicated for even for a fast magnetohydrodynamic shock. I'm not going to write them down in detail, but a version we use is from Wheatley et al in 2005. When you substitute this back into the characteristic equations you get a set of equations that describe M as a function of A of X if A of X is known. So if A of X is known and quantities ahead of the shock are known, B is the magnetic field, F would be the angle of the magnetic field to the oncoming shock, which is one way of specifying a parameter, then the system is closed and can be integrated numerically. You can't do much with this analytically because the shock jump conditions are more complicated for MHD. For one dimensional cylindrical convergence, A is proportional to R. So you know A is a function of distance in advance. You can integrate these numerically and get a solution for converging MHD shock. So just as a validation of the geometrical shock dynamics for MHD, we did a comparison of the GSD, geometrical shock dynamics, with computational magnetohydrodynamics where you solve the full MHD equations for an ideal oil or gas using a seven-wave method of some tanny. This is just an artificial example problem over a shock. Cylindrical shock converging in the presence of a spiral magnetic field form of a line current plus a radial magnetic field. Now, of course, the radial magnetic field can't exist because there are no magnetic monopoles, but you can put it on a computer and you can calculate the solution. And then GSD, MHD comparison, these are four different properties. The dotted lines are the computational MHD. This is done from full ideal MHD. The sort of lines are the GSD computations, and they seem to agree with each other reasonably well. So GSD does a good rendition of this problem. If you want to go to two space dimensions, then you don't know A as a function of distance in advance. So you have to go to a kinematic description in which the shock moves along its normal and one constructs as the shock moves a system of rays which are the orthogonal coordinates to the shock itself. And Wittem wrote down a set of kinematic relationships for the relationship between a surface moving normal to itself with a given Mach number and the ray coordinates itself. So this is purely kinematics. To close this, you need the area as a function of Mach number, the magnetic field angle, and possibly the distance if you know this. And this will give you a closed system that you can solve numerically. I'll write down a version of this a little later. We solved this using a numerical method due to Schwenderman and adapted this to the present MHD formulation. The physical problem we're interested in, which is the second part of my talk, is the evolution of a plain perturbed shock. This gives one example, both for a gastronomic case and an MHD case, the shock is moving from left to right. It's initialized. The shock is initially flat. Boundary conditions are periodic in the e to direction, and the shock is perturbed at t equals zero by just a perturbation of Mach number along the shock, which is periodic. And the shock has a Mach number of 1.8 in this example. And the perturbation is actually quite large. The solid lines show the shock shape as the profile, shock profile moves from left to right. Initially, and it's amplified just so you can see it by five folds, just for display. The shock is initially flat. Everything is smooth. But after a finite time, what is called a shock-shock develops, which is a kink in the shock wave and a discontinuity in the shock shape. And also discontinuity in the Mach number along the shock. The reason this develops is because waves propagate along the shock. These reinforce, just as they would in a Berger's-like equation, and they themselves form a shock along the shock. This was referred to in the parlance of shock dynamics as a shock-shock. So we can see these simulations show that from smooth initial conditions, a kink develops in a finite time, at least for this example. And you can see these contours, the contours of Mach number. So this contour gives you the Mach number at the position, this position when the shock passed that position. And you can see, you get this triangular cellular pattern. MHD shock behaves rather similarly, except the symmetry is broken by the presence of the magnetic field. But you can still see the presence of a shock-shock in a finite time. So we conducted a series of numerical simulations for both gas dynamics and MHD shocks. And we're interested in the time to formation of a shock-shock. And this is shown on this axis versus asylum. And we did a number of different cases, both for gas dynamic shocks and also for MHD shocks with varying magnetic fields and varying angles of the upstream magnetic fields. Here are some results. This is a empirical suggestion that, it's a line that suggests that the time the shock formation goes inversely is proportional into the amplitude. We can't get to very tiny amplitudes, but it suggests that no matter how small the amplitude, eventually a shock-shock will develop. So shock-shocks will not remain plain. This is not a stability issue. The amplitude, as you can see earlier, doesn't really grow, but it develops a kink which persists. So this is an issue of nonlinear dynamics. You cannot capture this with a linear approximation. So this also exists for cylindrically converging shocks. Here is an example of a sector. This is done in a sector of a cylindric, the shock is initially a circle, but it's perturbed along its arc length by a Mach number perturbation, which is periodic around the circle. The left shows the gas dynamic shock. The right shows an equivalent fast magneto-hydrodynamic shock, which is collapsing onto a line current at the origin, which it produces azimuthal field lines. These look quite similar, but they are dramatically different in one important respect, and that is that in the gas dynamic case, these contours are coloured by pressure, and the gas dynamic case converges to a singular pressure at the origin, which is the reason, or one of the reasons for the inertial confinement fusion, creates a pressure-temperature hotspot. But in the case of the MHT shock, this doesn't happen. The pressure actually goes back to the ambient pressure, so you do not get a hotspot, and that's a problem for using, if you wanted to do an inertial confinement fusion device in the presence of magnetic fields. The shocks behave very differently, and I just mentioned that in passing. So this is the shape, oh, you can see earlier also that you get the formation of shock-shocks. They're a little bit different, but not much different for the two cases, but the pressure behaviour is quite different. So you get the formation of a shock-shock in a finite time. This shows the gas dynamic shock implosion in the MHT, a perturb shock implosion on the right. Again, quite similar. This is the final stages of collapse, but you can see that there is a pressure singularity here, but no pressure singularity in the case of the MHT shock. So this formation of shock-shock seems to be, as far as we can see, a universal nonlinear phenomena occurring in these plain, perturbed shocks. So we wanted to see if there was an analytical basis for this. And so in order to do this, this is part three, we looked again at the equations of GSD for plainer shocks. This is the equation that I didn't write this down earlier, but this is the x is just the shock position, and this is just a kinematic description describing a shock which is moving in the vertical direction along its normal. This is just a kinematic description of a curve that moves along its normal with velocity A0 times m, where A0 is the speed of sound in front of the shock. This is just a kinematic relationship. You can write down a complex variable version of this where z is x plus Iy, and s is just a unit vector along the shock. n is a unit vector normal to the shock. You can write down a parameter along the shock, and t is the time parameter. This is just kinematics. You can express s in terms of the derivative of beta, and you can derive a nonlinear kinematic equation, which is the equivalent of the 2D GSD. It's not closed because you don't know q as a function of A, where q is going to come from the MA relationship, which is given by this expression. We wrote this down earlier for the shock moving down a tube of varying area A, and the area A can also be expressed in terms of b. This is just a kinematic equation. When you substitute into it the strong shock limit for which q of A, the Mach number as a function of area is A to the 1 on p, we mentioned this earlier, then you get a description of two-dimensional geometrical shock dynamics, and here I've written it down in terms of a shock that's moving in the vertical direction with constant speed, which is given by this, plus some perturbation. If z was zero, this would just be a plain shock motion, and it would just give you d z d tau. It's a constant, tau is just a rescaled time. This equation is a closed equation. It's a nonlinear wave equation, but it's highly nonlinear because of this fractional power of the derivatives of z and z star, where d is just a d d beta. We want to see if we could do an asymptotic expansion of this equation, and the approach we're going to use, this is just the same equation written down again, p is just a constant, but it's a known function of gamma. I should also mention in passing you can write down an equivalent equation, a similar but different equation for the weak shock limit. I'm just going to discuss the strong shock limit. So we're going to see if we can do anything analytical with this. And the starting point, this is a fractional power, it's a binomial series. Now, this is a complicated piece of analysis. And I'm just going to go through, sketch it in without going into too much detail. So the starting point to this equation is to do a binomial expansion of this guy, in terms of a quantity p, which I'll describe in a moment, p is just given by this, just a series of binomial, this is just a binomial character, we want to then use a Fourier representation. And what we're interested in is the behavior of these Fourier coefficients, uniformly as a function of time, when m is large, you can form the various Fourier convolution sums. And this analysis gets very complicated. But fortunately, with the use of symbolic manipulative tools like mathematics, you can simplify the analysis. Substituting this back, doing an expansion to order k is equal to two gives a sequence of ODE's for the Fourier coefficients. These are quite complicated, but you can solve them using a recursive method. One needs an initial condition given by this expression where b, the b1 is the definition of essentially the amplitude of the initial perturbation. And the asymptotic solution of the ODE's can be carried forward. You need a number of layers to do this calculation, I'll just mention them very briefly. One starts with an unzut for the Zm, this just says that the Fourier coefficients are going to have exponential convergence for at least a finite time. From this, you can derive an expression ODE's for the Rms. It's rather complicated, I'm just going to sketch it in. But these can be, this gives a sequence of ODE's for the Fourier coefficients, which you can solve recursively by which I mean the right hand side turns out to be just a function of the coefficients, the Zr's for m, for smaller values of m. And you can solve these recursively up to a set of amplitude coefficients, which I have been written down. And you might get stuck there, except there's a very clever generating function method developed by Derek Moore in 1979 to look at singularity formation on vortex sheets, which is directly applicable here. After one does some analysis, you actually get a solution which gives you the Fourier coefficients for large m uniformly as a function of time, it's a very complicated expression, but all these constants you see here kb, c plus q are known as functions of gamma. So you actually get a hard result. Note you get a m to the minus five half roll off multiplied by a exponential to the power m. So this Fourier series retains exponential convergence provided the argument of the exponential is always negative. In other words, Fourier series loses analyticity when this quantity is negative. And so this gives you the critical time for which the Fourier series is valid. At this critical time, you get an m to the minus five half roll off, so that the Fourier series is no longer analytic. Loss of analyticity in a finite time, inversely proportional to the amplitude with all constant snows. I'm sorry? Yes, it turns out to be a maximum. Yes, you can show it's positive. Yes, that's correct. You've got to be careful. There's a lot of detail that I just can't go into in this, but you have to be very careful about all of that. So you get one on epsilon. In his study of vortex sheet dynamics, Derrick Moore found that the behavior of the critical time was log epsilon, not one on epsilon. The difference here is that he had an extra e to the minus m in his expression, which led to a log. We have oscillatory expressions in place of that. So we end up with a one on epsilon instead of a log epsilon. So this result tells you that the motion of the shock will lose analyticity in a finite time given by a proportion of one on epsilon. We claim this also holds for MHD flows, even though this analysis has only been for gas dynamics, is also these are the shock jump conditions, at least for a magnetic field, which is parallel to the unperturbed shock. This leads to just an order of silo over m squared correction. So we believe this holds for at least strong MHD shocks. Comparison with our earlier results, here our earlier numerical results. Here's our analytical result with no fudge factors. We underestimate the time for shock-shock formation. And the reason for this is shock-shock formation is actually a derivative discontinuity. M to the minus five-half roll-off gives a curvature singularity and not a discontinuity in Mach number. So we think that loss of analyticity is a precursor to the formation of shock-shocks. And it's not actually the formation of shock-shocks. And we think that's probably the reason for the consistent underestimate. How am I going for time? I've got five minutes left. So the analysis gives us a hard result, tells us that no matter how small the amplitude of shocks, your initial perturbation is that eventually the shock will evolve to loss of analyticity. And then we believe this is a precursor to the formation of a shock-shock. A shock-shock is of course physically the formation of a triple point on a shock. Okay, so what does this singularity look like? Here is a comparison between the Fourier analytical methods. Once you know the Fourier coefficients, you can sum them to get the actual shock shape. This is a Mach number distribution in mind of the perturbation Mach number distribution. This is tau increasing. This is just for one case. Fourier analytical GSD numerical. And they agree reasonably well. There are some differences, but they agree reasonably well. Here, this is just some finite time up to 25. I think the predicted critical time for this is 31.614. And here, what happens here is that the Mach number becomes almost discontinuous, but not quite discontinuous. And it's similar to the GSD result, but here the shock-shock hasn't quite formed. If you look at the actual profile, Mach number, the shock profile, the shock's moving up. And you can see here is a shock profile with tau increasing at an almost but not exactly develops a kick. So this is a shock profile. This is a graph which shows on the left the slope of the shock with tau increasing. And the last one right at the critical time shows essentially a vertical slope at the point where we believe the ameloticity is lost. And this one here shows the curvature. Although we can't technically prove there's a singularity in curvature, it looks like that the loss of ameloticity is because the shock curvature becomes infinite while the slope remains continuous, but at a small time after that the shock-shock develops. So this is all I have to say. What we've done is an extension of GSD to MHD shocks. We've got some numerical solutions of the MHD or the GSD shock for both planar and converging shocks. For 2D perturbed shocks, the perturbations and the shock shapes are initially all smooth, but in a finite time a shock-shock develops. It appears to be inversely proportional to the initial amplitude. We then developed a GSD formulation in complex coordinates. And it could be, I gave an example with strong shocks, which you can do an equivalent one for weak shocks. You can actually do the general shock case also, but you have to do it more numerically. You can't write it down as a closed equation the way I have done, although it can be done implicitly. We've developed an analytical solution for perturbed planar shocks using a Fourier representation. This gives an asymptotic analysis for the Fourier coefficients for large M-uniformly in time. It gives to a time to singularity formation to leading order, which is inversely proportional to amplitude. We think the analysis also applies at least to strong MHD fast shocks. The loss of an analyticity corresponds to a curvature singularity, and we note that Tc, underestimate, I've written Tb, it should be Ts here. We think loss of an analyticity is a precursor to shock-shock formation, which physically would be triple point formation. The bottom line is perturbed planar gastronomic shocks do not remain plain. They develop according to this theory kinks in a finite time, no matter how small the amplitude. This has all been based on approximation. Shock dynamics, to verify we need to do something a bit stronger, and the current work is quantification using a shock-fitted Euler code. That's all I have. Thank you.