 So, we have a short discussion now, so I'm kind of combining some new work today with some very classical ideas of yesterday we talked about very simple source model. Today I want, now in the last half hour I want to go in a slightly different direction and talk about where I see a lot of opportunity for the future. This is going to be work that was done with Kyle Anderson as part of his PhD thesis and also now with Inche who's here with us today. And the idea being that volcano geodesy alone can only tell us a limited number of things about what's happening in the subsurface. Talked about yesterday the fact that we can get something about the depth of the magma chamber, some crude ideas about its shape. We can get the volume increment but not the volume. So there are things that we learn that are very useful and we can see inflation preceding many eruptions. But we don't get all the things that we would like to do and we're not really in a position to use them in a strictly forecast sense. And so what we have been working on and I think there's a lot of opportunity for the future, sort of where our group is going to put a fair amount of emphasis, is to combine the deformation measurements that we make from GPS and INSAR with the kind of conduit models that Michael talked about yesterday. So I think this is hopefully going to be nice for you because it's going to merge some of the things that you've heard about. So again, just to emphasize, we've seen this before. For some reason I have to look this way, I can't help it. But we can measure deformation. We can measure these inflation-deflation cycles and now we see these exponential-like recoveries, 1 minus e to the minus t, type behavior in these cycles. So this is something that over the last few decades as GPS has become more mature, INSAR has become more mature, we're getting better and better at doing this in multiple places. At the same time, we've learned a lot about how physical properties change as pressure decreases as parcels of magma move towards the surface to all the things that Michael talked about yesterday. So the goal that Kyle and I set out to do in his thesis was to try to merge these two disciplines and apply it specifically to the dome-forming eruption at Mount St. Helens that Michael also addressed. So again, remember, we had an explosive eruption in 1980, a big-plenty interruption column, and then some two and a half decades later, this very degassed, day-side spine was extruded. So very similar kinds of chemistry, but very different type of eruptive behavior. GPS that we had at the time, this is from the USGS, Mike Lasowski, was pretty good. You can see as the eruption occurred that there was deflation, the GPS stations moved towards the summit of the volcano as we'd expect. We only had one continuously recording GPS receiver in place prior to the eruption onset that's Johnston Lynch Observatory. And you can see there actually was very little inflation going on, very rapid deflation early in the eruption. And then interestingly, relative to what I was just talking about, as soon as the eruption lasted about three and a quarter years, and very rapidly began re-inflating. So now I have this question, does that really mean there was recharge coming in? Or maybe this was this viscoelastic effect? I don't know yet the answer to that question. I'm going to show you some model calculations. Once the eruption was underway, more continuous GPS sites were installed, these are their time histories. So we have the temporal variation and deformation at a number of points. And we also have estimates of the extruded volume. This came just from taking digital elevation models done from photogrammetry and other methods and differencing them to estimate the volume of the dome that was extruded. And I'm sure almost everyone has seen this sequence, but just in case you haven't, I'll walk through. This began in the autumn of 2014. This is Mount St. Helens crater. The lateral blast went off to the north. This is the dome that was left over from the 1980s eruptions. And then this welt came up underneath the glacier and then pushed out against the south wall of the crater. You can see it pushing the glacier to the side over here. So from that, we get the volume as a function of time of the extruded volume. So what, as I said, Kyle, Michi, and I now are trying to do is to build a relatively simple model of an erupting volcano, couple that to an elastically deformable medium so that we can compare to surface GPS measurements. And so there are going to be some of the features that Michael talked about yesterday. We're not going to get into some of the details at the scale of the bubbles. We're going to assume this is a relatively slow eruption. So we're going to assume that all these processes are at equilibrium. So we don't have to worry about the diffusion, the time scales of diffusion of volatiles into the bubbles, for example. But we'll have a magma chamber here from fitting GPS measurements. It's been fairly well established that it's more like a proletalypsoid than a sphere. There may or may not be recharge coming in to the bottom in some way that we're going to assume also is similar to the way I had it in the previous analysis. Magma moves up a roughly cylindrical conduit. And at some point, it completely freezes into this essentially solid day-side spine of extrudes up along a frictional surface, a ring fault system to the surface. And so things happen as the magma moves up. Exolved water and CO2 come out of solution and form bubbles. Some of that gas may escape out the sides or up the top. Crystals were formed. And as those crystals make up a larger percentage of the magma, it will eventually turn into something that's semi-solid and then slipping along these frictional interfaces. And we know that because there were earthquakes that were recorded that are interpreted as being occurring on these boundaries. So the idea is to see, can we learn more about the system from making a more physical model of what's going on? And remember, I always like to emphasize, our job as volcano geodesists is not to fit data. Our job is to learn about the magnetic system. So we can always fit the data. We could fit data with models that we know are wrong. You say, I fit the data. You say, so what? What have you learned? So the idea is we have to learn something. Now, of course, there's a trade-off. I can make models that become more and more sophisticated and have more and more unknown parameters. And if I get enough unknown parameters, I can fit anything. So there's a balancing act here. We want to make this as realistic as we need to, but not too complicated. So you have to watch the time, Ellie, and maybe start waving at me when I've used up too much time. So I think the details, in some sense, are not the critical thing here. It's the approach that I hope you take away. So this is from Kyle's thesis or two papers published. So we have essentially the same kind of equation for recharge as I had before. The time rate of change of pressure in the magnet chamber depends on the flux in versus the flux out. So we have magma coming into the bottom. We have flux that's erupting out at the top. We have our conservation of mass and momentum on a radially-average 1D profile. So this is exactly the equation that Michael derived for you yesterday. This is the continuity equation here, conservation of mass, one-dimensional conservation of momentum. Magma viscosity, as he emphasized, is a complicated function of dissolved water content and crystal content, the higher the crystal content, the much higher the viscosity. And then once this turns into a solid plug, we have friction on the boundary. And then we have mixture properties, say, for the density of the mixture of these different phases, a liquid phase, a gas phase, and a solid phase. We have to worry about things like solubility of CO2 and water. The reason we include CO2 is not because CO2 is volumetrically that important, but there are two reasons, one of which is that there could be measurements of CO2 discharged. That could be independent data that we can feed in. But also because when you have a single water phase, you get a step change, you have a discontinuity in the physical properties and the water first appears as an exalt phase. So not too important other than to say, as pressure decreases, the CO2 comes out of solution much more rapidly than the water, as was emphasized yesterday. Now, one of the things that Kyle did in his work was to build a relatively simple forward model of the conduit in the erupting system, compare that to observations, the observations being both the GPS data and the erupted dome volume, but then to do the inverse problem through exploiting something called Bayes theorem. How many people have seen Bayes theorem in a geophysics class? Oh, many of you, so that's very good. So what we have here is the so-called posterior probability, that's the probability of the model parameters, whatever our unknown parameters are, conditional on the data being proportional to what we call the likelihood, that is the likelihood of observing data D given model parameters and a prior probability. The prior is all the information we know about the system before we collect the data. And we can sample from this posterior distribution using so-called MCMC methods that I'll sort of elucidate on the nice slide that Kyle made. So we take an initial guess at the model, time step through the governing equations, analogous to the code you ran for your exercise yesterday. Of course, this is not an explosive eruption so we don't have fragmentation. This is a slow eruption where inertial effects are negligible. But you time step through until the eruption's over, then you make predictions, you calculate the observables at all the time steps, this would be the displacement of points at the surface, how much volume has come out, you can pre-compute Green's function is to do that calculation, compare that to the data. So you have predictions, you compare them to the data to get residuals and then you sample through that to something called the Metropolis Hastings rule to determine whether or not you retain a sample and then you just do this many, many, many times. So in order to do this, millions of times that forward calculation has to be something that you can run very quickly. The nice thing is you build up distributions then of the parameters that we're interested in. And so again, I don't wanna emphasize the details but I'd like to take away the general message here and this is a nice figure that Kyle made showing distributions of the parameters from what we would call the conventional volcano geodesy approach that is fitting ellipsoidal magma chambers or mogie sources or what have you, but we typically get, I'd emphasize this a number of times, you get some information about the depth of the magma chamber called at the centroid depth. So here's a probability density function of the depth. You get something about the shape, that's here the aspect ratio that is how much longer it is in the vertical direction than the lateral, right? And that came from the ratio of horizontal to vertical. The stock like, the more elongate is more effective at producing horizontal displacements. And you can get the volume change or the product of the volume or the total volume and the pressure change. Remember in the mogie solution we got P a cube. You don't get P separately and you don't get a cube separately, you get the product. So this is what you get from the sort of conventional analysis. When you build a more physics-based model there's a lot more going on. There's dissolved water, how much water is they able to depend on how much of that water turns into an exhal phase and it becomes a bubble that affects the compressibility of the system, compressibility determines how much volume change you have in the magma chamber versus how much is coming out of the top, okay? And it allows you to also place more physical constraints on some of the parameters. In other words, you all know this, I think we may have talked about it yesterday, there's this product of P a cube to occur but remember we said if the radius got too small the pressures would be unphysically large to explain the data. So we can start building in bounds based on physics and when you do that all of a sudden now look at all the different things that you have some hope of resolving and getting probability density functions for. So the same things we had before centroid depth aspect ratio volume change but now because we have these other physical constraints we can say things like what is the total dissolved water content? Because that affects the compressibility of the system is the recharge coming in from below. What is the total volume? So things that we care about as volcanologists were able to put better constraints on. Is the red line a product of P a cube? The red line I think is the maximum probability, maximum posterior probability. This should be, he had uniform bounds and I believe the box shows what the prior bounds were. Sorry? You're using net GVS expression for the thing that's going to put it in the space? Yes. Using the time dependent, full time dependent. You're all pretty big water and you need to keep that in the cold water or are you pretty involved in the same? Oh yeah, yeah, yeah, yeah. It's this model. It does have time dependent magma chamber pressure as magma's flowing out, the pressure's decreasing within the magma chamber. So it is a dynamic model. Okay, and then he was able to build posterior PDFs on something about the shape of the magma chamber. Okay, so now moving on to what we've been working on with Ichi lately, very influenced by a paper by Schneider-Rempel and Cashman that put in a little bit more physics than Kyle had in his analysis. In particular, we put in two things. We put in equilibrium crystallization, we put in gas loss and we put in a more natural rheological transition from the liquid to the solid. So in something like 10 minutes, I'm gonna go through this very quickly. And please save me or Yingqi if you want more details. So again, crystals will form and we're gonna take that, the solid volume fraction to be a function of pressure only, assuming equilibrium crystallization. So we're ignoring kinetic effects. Bubbles will form and gas can escape either vertically or laterally. They will have a transition from viscous flow to plug flow. We use the Schneider et al analysis. They calculated the solid volume fraction, the crystal mass fraction as a function of pressure using the thermodynamic program melts. So for a given temperature, we can calculate what fraction will be solid. And we're gonna assume that this is isothermal, that the radius of the system is large enough that we'll ignore heat loss for the most part. So as crystals form, the viscosity is gonna increase by a tremendous amount. We assume that we know that rheological law and that the total velocity, remember, we're doing a radial average velocity is the sum of a viscous part in a frictional part. The viscous part will look just like we spoke about yesterday where there's a shear stress that depends on the difference between the pressure gradient and the magnetostatic gradient. R is the radius of the conduit, inversely proportional to the viscosity. And then there's a rate dependent friction that could occur on the boundary which is proportional to the normal stress. And so what happens is that at great depth, the normal stress is high, the crystal content is low. And so all the flow, essentially all the flow will be viscous. As you get to shallower depth, the crystal content will get higher. The viscosity will get exponentially higher. And yet the normal stress is decreasing. So we'll transition rather naturally to the frictional flow of dominating. It's kind of a busy slide. We have continuity of both the liquid in the solid phase and continuity of the volatile phases. Both water is shown here in CO2. So this is the time dependent term. We're actually gonna look at only steady state solutions for the moment. So water can be either dissolved in the liquid or exalt in bubbles. And then there can be both vertical permeable flow once there's a connected network of pores or water gas can flow out laterally through the sides. And we're gonna assume that there's a percolation threshold that with roughly spherical bubbles that you have to get to a critical density, critical volume fraction of the gas phase before a flow can occur. So there are a bunch of parameters that we don't know. You could describe this in different ways and to determine what we think we know a priori and what we don't. This is just an example. First of all, show a steady state solution that Ying-Chi computed. And let's just go through a couple of these to show you the behavior. So here this is depth along the conduit and these are the volume fraction of the different phases. Yellow is the gas phase, blue is the melt, and red is the solid. Remember we're assuming equilibrium crystallization that depends only on pressure because temperature is fixed. So as the pressure decreases, we get to a point in the solid fraction really radically increases. There's a threshold when permeability becomes active and at that point the gas volume fraction starts to decrease because we can lose gas to the outside. And eventually the melt disappears and it's either solid or pores. So we transition from a viscous parabolic profile to plug flow where the upward velocity transitions from viscous flow here goes to zero and all the sliding comes along the frictional interface that occurs at a peak in the shear stress. This is the pressure profile here. So this kink in the gradient is where the shear stress becomes maximum. And so Yin-Chi computed this for different temperatures. So remember we're assuming temperature is fixed but of course the properties are functions of temperature. So if we were to decrease the temperature from 850 to 800, so the velocity, the solid volume fraction would increase more that would depress this transition to the solid plug to greater depth and slow down the flow radically because the viscosity would be much higher. So in the last five minutes I'll show you what Yin-Chi has been working on and that is to take other kinds of observables. So we're gonna look at kind of a steady state part of the eruptive cycle after the initial period where the eruption began. It's never truly steady state but we're gonna approximate it. The point we wanna make here is there are other kinds of observables that we can add to this problem. We've added a bunch of parameters that we don't know. We don't know all the properties of the permeability. We don't know some of the other properties, the frictional properties and so forth. But we wanna show you that we have more observables. We do know the exit velocity. We can do that from photographic measurements but also the USGS did this wonderful thing. They slung our GPS spider on top of the dome. And so we know the extrusion rate from those measurements. You can see the spider just moving along the dome from that we get the linear extrusion rate. We also know the volume flux from the differencing the digital elevation models. We also know something about the porosity of the eruptive materials, right? And that puts a constraint on permeability because if we start with a given amount of gas and we end up with five or 10% porosity then we know so much of it must have escaped to get out and leave us with this kind of process. So that becomes another observable. We know from patrology, the crystallization depth, patrologists look at the ground mass and can tell us something about the pressure at which the last liquid phase crystallized. And we're gonna assume that we know something about this transition from viscous flow to plug flow from the drumbeat earthquakes that people at CVO have analyzed. So again, this is just sort of a candidate of what we could do with these kinds of information. So this is actually not even using the GPS data. We're not using the time-dependent GPS data. Just these observables that I just mentioned and Ying-Chi does the Monte Carlo sampling and ends up with these, the blue here are the prior distributions and the red are the posterior and ends up with some distribution for say the magma chamber pressure, the total water content, the permeability constant on a log scale, percolation threshold, friction coefficient and so forth. Point being that we can now put in more observables, constrain unknown parameters about the system and then work towards going back now to the time-dependent system where we can put the GPS, the time-evolving volume flux and be able to build up more and more constraints on the system. So there's some interesting observations that come out of this rather than focus on that. I think in my view, the future is how do we build systems that can incorporate even more kinds of data? We talked about gas emissions. We haven't included gravity. We're gonna hear more later in the week about gravity measurements. Other kinds of geophysical observables, geochemical gas measurements. These models give us the ability to rationalize these different kinds of measurements in a self-consistent formulation. And I'll just leave you with the question. Is it possible to use this in a forecasting way? This is something a number of us are grappling with. I'm gonna put this out as a question. I don't know the answer to this. But of course, these models are fully deterministic. We have conservation of mass, we have conservation of momentum. If we know all the material properties and we know the initial conditions, or we can calibrate them, then in theory, put this in, you know, in italics, in theory, we should be able to forecast forward in time once an eruption has started and maybe say something about how it will behave. And I'll just show you one figure I put together a number of years ago, just as sort of a thought experiment using the GPS data. And this is based on Kyle's analysis, not including any of the effects that you put in. But we have here the volume extruded as a function of time on this axis. We also included the GPS data, although I'm not showing that. So this is just volume flux as a function of time. And what I did in this calculation was take the data only up to here. Do the MCMC sampling to constrain the initial conditions, the initial volume, the initial pressure, the initial water content, unknown parameters and so forth. And then randomly sample from that distribution and say, now it's deterministic once I know the initial conditions, I can propagate forward in time. And this gives trajectories of the system that would go forward. This particular model was so simple-minded that it could do only one of two things. So it's not useful in a real world context. It's a more of a gedonkin experiment. This system could only do one of two things. It could either go into a steady state eruption, as shown here, where the flux into the bottom of the magnet chamber balance the flux out. Or if there was no influx, you'd get an exponentially decaying flux that would eventually, asymptotically end at infinite time. Because temperature isn't changing, it can't freeze, there's no solidification that can happen, it only has one of these two cases. But you can see that based on the data up to one and a half years, it sort of seems to favor the ending eruption rather than the steady state eruption. And you can do this in a probabilistic sense. And in fact, we weren't the first people to try this. Larry Mastin and the people at the Cascades Volcano Observatory were actually doing this in quasi-real time. They were using a simpler model than the one I just showed you, but they were doing the same kind of thing, seeing whether or not, based on the data that they have, that it would go into a steady state eruption that would go on for many decades, perhaps, or it would asymptotically end. So I do think there is the potential for building models of this sort. And it remains to be seen whether they can be made useful. But ultimately, again, the usefulness of models is whether they predict data that you don't have. You can always fit the data that we do have, but if we could actually predict something into the future and do true tests, forward tests, I think this would be incredibly useful. So I think this is a big challenge for the community. Don't know whether, you know, how successful we will be, but I see this as a very promising area moving forward. I think I should stop here and take any questions.