Interactions of Buoyant Upwellings with Subduction Induced Mantle Wedge Flow

One of the primary modes of thermal-chemical transport in subduction zones controlling the growth of crust and the evolution of ocean-atmosphere system through geologic time is buoyant upwellings or diapirs. Typical process models developed from surface data depict vertical ascent paths in two-dimensional cross-sections through convergent margins. In this investigation we show that compositionally buoyant diapirs sourced from the slab-wedge interface have significantly more complex communication pathways as they interact with three-dimensional, timevarying mantle circulation driven by plates. Analogue fluid dynamical laboratory experiments, capable of representing the needed range in length scales (10 orders of magnitude), are used to characterize diapir dynamics for a range in plate-driven and buoyancy-driven parameters. Circulation patterns are recorded in traditional crosssectional and map-view planes using photogrammetry to measure deflections from vertical ascent and parallel to the strike of the trench imposed by the creeping wedge flow. Results show a) the style of plate subduction produces different diapir paths, with spatial-temporal complexity (95% non-vertical paths) b) up-slab and down-slab buoyancy fluxes lead to three basic surfacing modes in arcs and c) threedimensionality of subduction-driven wedge flow, specifically trench-parallel motion, leads to complex diapir-conduit interaction scenarios, with strong implications for geochemical models of arc volcanics . Diapir thermal models produce a complex range of evolution scenarios, from no-melt, to partial-melt and to full-melt outcomes depending on individual paths. ACKNOWLEDGMENTS I would first like to thank my advisor, Dr. Christopher Kincaid, for his unwavering support throughout my time at the University of Rhode Island and Graduate School of Oceanography. Without his interest in supporting students from diverse educational backgrounds I would never have even applied for a research based graduate degree program. His patient guidance led me through more than three great years of intense and unfettered learning for which I am forever grateful. His passion and excitement for exploring various problems at all scales using a myriad of complimentary techniques (the three legged stool) is truly inspiring for researchers and students alike. Thank you Chris, I really appreciate it. I would also like to thank my thesis committee members, Dr. Steven Carey and Dr. Dawn Cardace, whose suggestions and insight helped steer this project from inception and significantly increased the soundness and quality of the work. Additionally I would like to acknowledge my lab mates: Sara Szwaja, Christina Wertman, Kevin Rosa, and Loes van Dam, for countless discussions, suggestions, and sanity checks. Conversations with Rob Pockalny, Brian Savage, Aaron Hirsch and WHOI collaborators Mark Behn and Nan Zhang greatly influenced the outcome of this work. Finally I must express the utmost thanks and appreciation to my friends and family, especially Mackenzie for sticking with me and seeing this through to the end.

Tables 1, 2) rising into the most simple plate-flow scenario, or downdip-only slab motion. a) and e): Initial conditions for two rows of 6 diapirs distributed in a trenchparallel orientation within the velocity boundary layer (VBL) above the sinking plate.
The use of two rows roughly represents the buoyancy input proposed for different fluid-release mineral phase changes, a sub-arc and a rear-arc output (Ringwood, 1982). where 4 diapirs are distributed in a trench-parallel orientation along each slab segment.
As in prior cases, diapirs are initially placed in the VBL above the subducting plate.
Colored circles are used to construct distinct diapir pathlines for those initiating above the steep belt (b) and shallow belt (c). (b) Patterns clearly show that flow/diapir transport is from the deep wedge above the steep slab segment to the sub-arc region above the inner edge of the shallow slab segment. c) Without rollback, diapirs sourced from the shallower slab also return to the sub-arc of this plate through curved ascent paths (1. slightly up and away from the trench, 2. near-vertical through midwedge, 3. strongly deflected sub-horizontal and toward the trench in shallow wedge).
Grey shaded circles represent starting points along diapir particle paths.

Introduction
Subduction is a primary mechanism of mass and heat transfer between Earth's surface and interior and a major driver of large-scale mantle flow patterns. Subduction also drives processes ranging from planetary resurfacing, the opening and closing of ocean basins, and the differentiation and accretion of continental crust. Tectonics and subduction play key roles in regulating the global carbon cycle and climate on geologic timescales. On shorter timescales, more relevant to humans, subduction is a large concern as it is intimately linked to destructive geologic hazards including large mega-thrust earthquakes, tsunamis, and explosive volcanism [Stern, 2002;Rosen, 2016]. Much of the global population is concentrated along the continental margins and consequently in close proximity to subduction zones.
We investigate the interactions that buoyant mantle upwellings have with subduction-induced circulation in the mantle wedge using analogue laboratory models.
Here mantle wedge is used to describe the region between the down-going slab and the overriding plate, which itself contains a sub-arc plate (closer to the trench) and a back-arc plate (farther from the trench). This study focuses on the origin of volcanic material in subduction zones, and building improved connections for how material moves from the deep slab-wedge interface (SWI) to beneath both sub-arc and back-arc plates. This work aims to expand on our understanding of where source material for volcanic arc magmas originates, the lengths and shapes of ascent paths, time of ascent, and thermal implications for different ascent paths.
One means of transporting subducted material from the SWI to the lithosphereasthenosphere boundary (LAB) at depth is buoyant upwelling or diapirism. Diapirism within the mantle has been proposed for years within the literature, for example [Berner et al., 1972;Whitehead and Luther, 1975;Oxburgh et al., 1978;Marsh, 1979;Olson and Singer, 1985;Olson and Nam, 1986;Weinberg and Podladchikov, 1995;Kelly et al., 1997;Manga, 1997;Hall and Kincaid, 2001;Gerya and Yuen, 2003;Hasenclever et al., 2011;Marschall and Schumacher, 2012;Miller and Behn, 2012], and is an elegant solution encompassing the episodic nature of surface volcanism consistent with the material transport mechanism. The diapir model is a parsimonious answer to the question of how to move relatively large volumes of material around in a complex and dynamic tectonic setting.
The physics of producing diapiric structures is relatively well understood and described analytically by the Rayleigh-Taylor instability and has been investigated experimentally and numerically [Rayleigh, 1882;Nettleton, 1934;Taylor, 1950].
Diapirism culminates with the detachment and subsequent ascent of a buoyant spheroidal body from an unstable fluid interface. Any layered system containing a fluid overlain by another, denser one is unstable. Small perturbations at the material interface between the two fluids result in imperfections and opportunities for the less dense underlying material to grow rapidly and rise through the overburden.
The nature of the density contrast between the two fluids may be thermal, as is the case for classical stove-top or Rayleigh-Benard convection; compositional as is the case for crustal salt diapirism; or a combination of the two such as a "lava lamp", which may be the most appropriate and accessible analogue available for demonstrating mantle convection processes. In natural subduction zones, the density contrast between buoyant diapir material and ambient mantle wedge asthenosphere will be both thermal and compositional in nature. In this study we investigate purely compositional density differences.  Hall and Kincaid (2001), is used to allow for diapir-diapir interactions and coalescence by incorporating rheological weakening of the matrix fluid. Our models make use of a subduction apparatus that is well tested, and includes capabilities for both downdip and rollback subduction styles for a complete range of fixed dip angles [Kincaid et al., 2013;MacDougall et al., 2014;Szwaja and Kincaid, 2015]. This apparatus (Figures 1, 2, and 3a) also has the ability to characterize 4-D wedge circulation and related diapir ascent pathways when the geometry of the subduction zone varies through along-trench changes in dip angle producing slab gaps that are fixed or can vary through time.
Our experiments are designed to explore the deflections of diapir paths caused by variations in diapir parameters (buoyancy, morphology, source spatial patterns) and plate-driven background flow regimes (subduction rate, style, geometry  The end walls are even farther from the mantle wedge ROI.

Ambient Working Fluid and Glucose Rheology
The creeping asthenosphere is modeled using highly viscous (96 wt% solids) glucose syrup as the primary working fluid following similar studies [Kincaid and Olson, 1987;Bellahsen et al., 2005;Piromallo et al., 2006;Kincaid et al., 2013;Druken et al., 2014;MacDougall et al., 2014;Szwaja and Kincaid, 2015]. The viscosity of glucose syrup is primarily a function of temperature and water content and follows an Arrhenius type function [Schellart, 2011]. Both of these parameters are held constant throughout the duration of a single experiment and between different individual experimental runs. At the small strain rates considered in these experiments any stress dependence of the glucose rheology is considered negligible.

Externally-Forced Kinematic Approach
Flow within the working fluid is induced kinematically through the use of a belt drive system (Figures 1b, 3a). Two continuous rubber belts representing infinitely strong subducting slab segments are driven by high torque variable speed DC electric gear motors. Flow in the wedge is achieved by viscous coupling of the working fluid to these belts, facilitating the development of a thin viscous velocity boundary layer (VBL) at the SWI. This system has been shown to reproduce mantle wedge flow fields that have been observed in fully dynamic (internally forced) analogue models [Kincaid and Olson, 1987;Richards and Griffiths, 1988;Kincaid and Griffiths, 2003;Funiciello et al., 2006;Kincaid et al., 2013;Druken et al., 2014]. The benefit of kinematically prescribing plate rates is an increase in experimental reproducibility, minimizing the variability between identical runs that internally forced dynamic analogue models suffer from.
Given that plate rates are the most well-constrained parameter of the subduction system in comparison with rheology, composition, and structure, an externally-forced kinematic modeling approach is justified. A maximum bound for the amount of energy put into driving this system can be derived from the maximum electrical draw of the drive motors, which is about 4.6x10 3 Joules for the duration of a typical experiment.

Degrees of Freedom
The custom-built subduction apparatus ( Figure 1) allows for control of two independent slab segments and their respective convergence rates (downdip motion is referenced by the variable Ud). Because the two belts are controlled independently there is an added capability of modeling along-strike dip changes, either fixed or varying over time (i.e. slab steepening or the opening of slab gaps). Deep belt rollers at the equivalent of the upper-lower mantle transition zone (~400 km depth) migrate forward or backward relative to a stationary roller (representing a fixed trench) to produce either slab shallowing or steepening. There is also control for variable trench migration or rollback speed and position. This becomes important for producing a toroidal flow component around plate edges (Figures 8, 10), significant for trenchparallel material transport along-strike. The upper surface boundary of the overriding plate can be left free, fixed relative to the trench using a small rigid sheet of acrylic, or prescribed kinematically using driven rolls of mylar sheeting to simulate back-arc extension. Experimental parameters for individual experiments are listed in Table 1.

Diapir Materials and Properties
A variety of materials were tested to model "cold" compositionally buoyant diapirs, see Table 1 for listing of experimental parameters and Table 2 for summary of measured experimental observations. The three distinct sets of experiments discussed above utilize three types of diapir material. The first experiments utilized arrays of solid rigid spheres. Buoyancy of individual diapirs was varied by utilizing various hollow and solid plastics and polystyrene foams (100-500 kg m -3 ) in a range of sizes (0.5-1 cm (15-30 km) radius). The primary motivation for this phase was to more fully map flow, especially near plate edges and the boundaries of the ROI, with variable mantle wedge geometries and convergence rates. This phase was also informative in assessing which regions would be conducive to successful diapir formation and detachment in subsequent fluid injection experiments and determining the optimal location for the injection point source. Solid rigid spheres were ideally suited for this phase, minimizing experimental waste before moving to more dynamically-realistic and experimentally-challenging fluid diapirs. These rigid spheres also eliminate the effect of diapir deformation or deviation from spherical shape that malleable fluid diapirs experience.
Regularly spaced arrays of the rigid spheres were embedded in the viscous boundary layer before entering the trench and subducted to diapir initiation and detachment depth (Figure 2a) [Miller and Behn, 2012]. Slab motion was then briefly paused to allow for a small amount of buoyant ascent and detachment of the rigid spheres from the viscous boundary layer at the SWI to initiate before subduction and instantaneous mantle wedge response was continued. Without the brief pause in plate motion rigid spheres remain embedded in the viscous boundary layer and fail to detach from the SWI and are not successful in ascending to the surface. This mapping phase allowed for the examination of the effects of initial position on diapir trajectory and surfacing location as well as determination of where probable diapir initiation and detachment may occur within the tank, as the potential depth range for diapir initiation (30-70 km) and detachment (90-170 km) reported in the literature are not wellconstrained [Hall and Kincaid, 2001;Gerya and Yuen, 2003;Miller and Behn, 2012].
In Air was also used, but is not reported here, see Appendix 1. Air in the lab has a density of approximately 1.180 kg m -3 given a temperature of 23 °C, 25 m elevation above sea level, and 70% relative humidity. The glucose working fluid has a density of 1420 kg m -3 , yielding a density contrast of -1418.82 kg m -3 , or nearly 100%. This large density contrast was explored as an end member most likely to produce ideal vertical ascent. The calculated ideal terminal Stokes velocity of a spherical air diapir with a 5 mm radius in quiescent fluid is about 3 cm min -1 (Figure 4).
It is essential that diapir fluid is not able to contaminate the large volume of working glucose fluid, as the changing of such a large volume of highly viscous fluid is a major undertaking. A glucose vacuum system was developed for removing diapir fluid that utilizes a heated wand to precisely remove the lower viscosity material. A procedure was developed whereby all weakened material was removed between experimental runs and the tank and wedge ROI were allowed to homogenize and equilibrate.
Diapirs are tracked throughout an experiment, from introduction until they either successfully surface through the wedge or become entrained deep into the working fluid exiting the ROI. A diapir is deemed successful if it forms, detaches, and transits the entirety of the mantle wedge from the slab-wedge interface to the lithosphere-asthenosphere boundary. Unsuccessful diapirs are therefore those that fail to form, detach, or overcome the advection of the background wedge flow and become entrained in it. We do not report on the fate of unsuccessful diapirs here unless otherwise noted for exceptional results.

Scaling
Scaling between the analogue experiments and the mantle of the Earth is achieved through the dynamical similarity of the Péclet and Prandtl numbers of the two systems [Kincaid and Olson, 1987;Kincaid and Sacks, 1997;Kincaid and Griffiths, 2003;Kincaid et al., 2013;Druken et al., 2014;MacDougall et al., 2014;Szwaja and Kincaid, 2015]. The Péclet number is a dimensionless ratio defined by the relative contributions of advective to diffusive heat transport (eq. 1). The Péclet number of both the laboratory and natural systems is on the order of ~15, suggesting advection is slightly more important than diffusive heat transport, but the scales of the two competing transport mechanisms are not disparate enough to warrant eliminating either term. The Prandtl number is defined as the ratio of kinematic viscosity to thermal diffusivity (eq 2.), which is very large for both systems, often being approximated as infinite for the mantle. Constant thermal diffusivity is assumed for both systems, 1x10 -2 and 1x10 -3 [cm 2 s -1 ] for mantle and lab respectively, and all of these experiments are carried out in an isothermal environment with no external heat added. We choose a characteristic length scale of 750 km to equal the total slab width of 25 cm yielding a length scaling factor of 30 km cm -1 , a time scale factor 0.5 min Myr -1 , and velocity scale factor of 1 cm min -1 to 1.5 cm yr -1 .

Cameras and Image Processing
Time-lapse photography is used to obtain flow field and diapir displacement data via two orthogonally oriented Ultra High Definition (4288 x 2848 pixels) DSLR cameras (Nikon D90). One camera was positioned at the southern wall facing northward through the wedge in a typical vertical cross sectional view with the slab dipping to the east. The second camera was oriented above the center of the trench looking down in map view (Figure 1a). Both cameras were fixed with the trench always yielding a fixed trench-centric reference frame. These cameras are interfaced through a computer for control and timing, and raw images are immediately archived to disk for later processing and analysis. We shot at a 5 or 10 second frame interval depending on the specific experimental setup.
Basic image processing techniques are used to extract relative motion between frames and diapir transport paths. A custom image processing pipeline was constructed using freely available open source software, namely Python, the SciPy/NumPy stack, and OpenCV. The pipeline consisted of a linear routine of basic operations (color space transform, histogram equalization, background model estimation, blurring, thresholding, segmentation, feature identification, and tracking).
Within the image processing pipeline we also measure other feature descriptors such as diapir radius, eccentricity, major and minor axis azimuth, rotation, Hu moments (scale, rotation, and translation invariants). These additional parameters help in identifying and keeping track of specific diapirs and dealing with occlusion, as well as filtering out misidentified parcels and other noise. This pipeline can be run in batch on all images from a given experiment producing 4D (x,y,z,t) Lagrangian position vectors for diapirs in each frame. Background mantle wedge flow field velocity structure is obtained using common FFT peak shift and interrogation window based DIC/PIV techniques [Thielicke et al., 2010].

Results
These experiments represent the complex interactions of two distinct mantle processes, background plate-driven circulation and the buoyant ascent of diapir material. For the latter, the ideal Stokes terminal ascent velocity is calculated for each of the diapir analogue materials used (Figure 4). Three basic regimes are relevant to these results in cross-section, 1) rapid downward flow of wedge material directly above the slab (we refer to this as the velocity boundary layer (VBL)), 2) the relatively sluggish interior of the wedge, far from the upper and lower boundaries, and 3) the shallow return flow of wedge material towards the trench and the upper and extreme inner, or apex region of the wedge, (Figures 6, 7). Flow regime 3 is important because it controls how material is brought to and along the underside of the arc/back-arc plates, exerting an influence on thermal-chemical exchange processes from the mantle to the crust.
Prior 4-D subduction laboratory models [Kincaid et al., 2013;MacDougall et al., 2014;Szwaja and Kincaid, 2015] show that downdip and rollback subduction produce very different shallow return flow patterns. Figure 9 shows  Figure 29). This pattern of steep deep-to-shallow dip ascent trajectories and linear surfacing patterns above the shallow dipping segment at greater trench-normal distances is enhanced when a slab window forms in rollback subduction. Figure 14 shows that toroidal flow around the retreating plate results in a much more confined trench-parallel transport of diapirs from steepdip wedge to shallow dip wedge (Figure 14b). Similarly, the toroidal flow around the retreating shallow slab results in a very pronounced linear distribution of surfacing diapirs, trending along a trench-normal direction (Figures 14c, 29). Interestingly, the surfacing times in this case do not reflect a simple, unidirectional age progression.
Diapir surfacing times reflect the more complex deeper interactions with flow resulting in more spatially heterogeneous surface arrival times to the base of the plate.

Source Processes and Diversity of Diapir Ascent Paths
A second set of experiments utilizes a continuous release of buoyant fluid into the VBL above the down-going slab to develop relationships between processes operating in this source region and the eventual ascent and surfacing patterns for diapirs. We begin simply, following the method of Hall and Kincaid ( 2001), focusing on release from a single narrow region above the center of subducting slab. This work extends that of Hall and Kincaid (2001) by testing the influence of both downdip and rollback modes of subduction (Figure 2b). Another simplifying factor is the use of a dry ethanol fluid for diapirs. Ethanol and other compositionally dry buoyant diapirs do not exert any rheological weakening on the wedge as they contain no water or heat to alter the matrix viscosity. There is no preferential fabric imprinted in the wedge matrix or weakened conduit features left behind by these dry diapirs that would allow two separate bodies to physically merge into a single morphology [Hall and Kincaid, 2001]. Ethanol therefore makes an excellent next step towards increasing model  Figure 16 highlights the arcuate ascent trajectory that is commonly observed, generated by a) rapid downward drag in the VBL away from the trench, b) slow, near-vertical rise into and through the midwedge and c) deflection towards a sub-horizontal trajectory during entrainment into return flow towards the trench and base of the sub-arc plate, with different mixes of trench-normal and trench-parallel motion depending on parameters such as subduction mode and geometry (Figures 6,7).
In these along-centerline release cases, the diapir source and surfacing locations are roughly aligned vertically because they remain within the same trenchnormal, cross-sectional plane (e.g., limited north-south offsets due to minimal trenchparallel flow). However, a number of key processes can occur at the base of the wedge within the VBL that displace surfacing locations either towards or away from the trench, relative to the deep diapir source region. Figure 17a shows  (Figure 18a-c). The rapid expansion and departure of a Mode II diapir can also occur due to up-slab flow of another diapir (Figure 15b,c). Here the sudden increase in size/buoyancy causes these diapirs to detatch earlier in its formation and surface closest to the trench. A third mode (III), where smaller groups or clusters of diapirs form down-slab to produce an increase in rise rate, is very repeatable in cases both with ethanol and rigid diapirs. Figure 19 illustrates the time evolution of a clump of smaller diapirs. Smaller diapirs have slower ascent rates, remain longer in the VBL and are dragged deeply into the wedge. They tend to separate and stall in the lower, mid-wedge where weaker buoyancy is balanced by weaker downward drag (Figure 19a). Once 2-3 small diapirs enter this region of diapir stagnation, a cluster forms with an integrated buoyancy that significantly exceeds downward drag (Figure 19b). These rapidly rising features, starting from deep in the wedge, produce near-vertical ascent trajectories and surfacing locations that are furthest out from the trench (~540 km) (Figures 17a, 19c).

Trench-parallel Flow and Diapir-Diapir Interactions
A third set of experiments explores how distributed along-trench buoyancy sources combine with trench-parallel return flow to produce complex diapir-diapir mergers and interactions. An along-trench distribution of diapirs is expected to occur naturally from instabilities arising from the sheet of buoyant material believed to exist above the down-going plate. The work builds from experiments by Hall and Kincaid (2001) that utilized a water-bearing ethanol solution (vodka) that produces rheological weakening in the wedge, diapir-diapir mixing and networked conduit features in the wedge. This is opposed to rigid spheres or dry ethanol diapirs that cluster, but do not coalesce into a single coherent body. The wet ethanol has a laboratory density contrast that is also closest to mantle values. A hydrous ethanol diapir with 1 cm radius has an ideal Stokes ascent rate in glucose of 4.4 cm/min, neglecting the effect of water on the matrix viscosity (Figure 4). Diapirs of this type tend to exhibit less trench-normal deflection than their dry counterparts because water at the diapir-mantle interface reduces the viscosity at that boundary, consequently reducing the drag by lubricating the surface.
The third experimental phase consisted of adding additional point sources along strike to better mimic a line source. The sources were configured either as a set of three or four sources centered on either the whole plate centerline (E_13, E_14, V_5, V_6) or on the southern plate segment (V_7-10). The spacing between point sources can be set to represent the instability wavelength, often thought to control the regularity of volcanic vent spacing (~70km) along strike observed at many subduction zones [Marsh, 1979]. Employing multiple point sources along strike is motivated by the fact that we are trying to model a 3D layer with thickness and structure becoming unstable with point sources. When a buoyantly unstable layer intersects the instability initiation zone it is transformed into a line source of buoyancy [Marsh, 1979]. The technique of approximating continuous line sources or sinks with points has been used to study other types of subsurface flows of importance to hydrology and petroleum engineering [Weijermars and van Harmelen, 2016].
Hydrated ethanol diapirs in these experiments reproduce some of the behaviors observed by Hall and Kincaid (2001)in that they do produce rheologic weakening in the wedge creating preferential conduit paths and networks for subsequent diapirs to follow. The weak zones are passive features with minimal buoyancy so they are completely entrained in the creeping matrix flow and redistributed throughout the wedge. This not only reinforces the hypothesis that preferential weak zones can exist and provide a means of rapid ascent, but shows that those conduits are not fixed in In some cases these swollen surfacing diapirs can begin deforming within the stronger return flow below the overriding plate to produce an extensive under-plating feature that tends to elongate into the direction of flow. Figure 24 shows such a progression for experiment V_9, which is similar to V_7 but with 40° dip angle, downdip only subduction. One difference is that this case has a no-slip (zero vertical alignment of buoyancy sources and rapid upward assimilation. An interesting process seen in experiment V_10 that has a no-slip condition beneath the back-arc plate, and therefore slightly enhanced vertical flow within the wedge (Figure 27c), is that drifting over-thickened tails remnant from prior upward diapir assimilation events, may capture smaller diapirs that are otherwise stalled deeper in the wedge, allowing them an access path to the surface.

Slab Steepening -Experiment D_14
During the course of this phase of experiments we ran into a few special and noteworthy cases that warrant reporting. The first of these, experiment D_14, used a complex subduction geometry with active slab steepening employed for the southern plate segment only. This could be viewed as a slab tearing with the northern piece continuing on the initial shallow path while the southern half drapes back into the mantle. This is not reflective of any particular subduction zone on Earth but rather an exercise in experiment design that the apparatus is able to produce.
Both plate segments initially dipped to the east at 37 degrees. The initial arrangement of diapirs was the same on each plate segment, one trench-parallel row of 4 diapirs. The diapirs were subducted and convergence briefly paused to allow for detachment from the viscous boundary layer. When subduction resumed, active steepening of the southern plate segment at a rate of ~1 degree of steepening for every 333 kyr was also started contemporaneously and ceased when the southern slab had The diapir from the southern group that originated closest the southern edge and hence farthest from the shallow plate segment had a much different trajectory (magenta circles, Figure 28). The other three diapirs from this side were advected nearly horizontally along strike toward the shallower fixed dip northern plate segment.
The southern-most diapir, however, was not advected close enough to the edge of the shallower plate before a bulge of material was extruded out of the newly formed gap between the two plate segments. This extrusion of material at depth pushed back against the lateral forces pulling in toward the shallower plate and caused this diapir to become nearly motionless at depth for a long time ( Figure 28).
The residence time of that specific diapir was 18.33 minutes (36.7 Myr) , nearly twice as long as the average of all of the other successful diapirs from that experiment ( Table 2). This kind of stochastic behavior is not predicted by conventional This corresponds to ~70% of the total trench width, or 525 km of lateral along strike transport.
This experiment produced an evenly spaced linear track of LAB surfacing locations at the centerline of the total slab perpendicular to the strike of the trench. Although the spacing between diapirs at the LAB was regular (~130 km), the source regions for subsequent surfacing diapirs were not characteristic of this surfacing order. That is to say material was arriving at the base of the lithosphere at regularly spaced intervals but was not diagnostic in predicting where a subsequent or preceding diapir would have been sourced from. A diapir originating from the slab edge could arrive at the LAB following one from the middle of the slab or vice versa. This has consequences for any petrogenetic relationship along a linear volcanic track because the vents are not sampling material from the same region of the slab-wedge interface, and any magmatic lineage would not be completely obvious.

Up-dip subduction channel flow -experiment V_3
In this experiment diapir material was observed migrating back up dip along the SWI, which is another proposed route for buoyant material to travel, especially if the SWI is weak [Gerya et al., 2002;Marschall and Schumacher, 2012;Harlow et al., 2016]. In this case a diapir burst at the surface of the working fluid and spilled into the trench (a behavior produced only by wet ethanol) before being evacuated with the vacuum. The hydrous ethanol was actually drawn down and subducted creating a thin film along the SWI partially lubricating the interface and partially decoupling the wedge fluid from the slab. When this film intercepted the injection source any fluid emanating from the source immediately preferentially flowed up this weak zone instead of contributing to the formation of a new diapir. This same effect could work in nature, but it has to be relatively localized, as flow in the wedge would cease if the slab and wedge asthenosphere were not viscously coupled.

Unsuccessful diapirs
Diapirs in the laboratory setting can fail to reach the LAB for a few different reasons. Unsuccessful diapirs in the solid sphere experiments were generally those that did not lift off or detach from the viscous boundary layer at the slab surface substantially during the brief pause in plate convergence. Inconsistencies in detachment are an unfortunate reality of this method, but mostly mitigated by the use of multiple diapir tracking targets and regularly spaced arrays. If one diapir is strongly embedded in the viscous boundary layer and fails to detach, the neighbors provide us with some of the same information. Not all material that enters a subduction zone is expected to successfully resurface as diapirs so this does not negate our conclusions.
Solid diapirs were also sometimes unsuccessful for experimental configurations that had a greater vertical velocity component in the wedge caused by steeper slabs or faster convergence rates. Since the rigid spheres have fixed buoyancy (they do not change size or density), the primary adjustments must be made to the subduction side of the force balance acting on a diapir. In the fluid diapir experiments we had the ability if needed to alter the buoyancy flux to create successful diapirs while maintaining a desired model configuration. We also experimented with populating the underside of the plate with solid spheres as well to see if they would be brought into the wedge through a slab gap or around the edge by toroidal rollback flow. None of these sub-slab diapirs were successful, so we do not expect any sub-slab material to enter or influence the wedge. In general for all the cases, diapirs that detach farther down dip have less of a chance of successful ascent. Unsuccessful fluid diapirs were those that were small relative to the average size for a given experiment. They failed to develop fully and detached from the buoyancy source prematurely before becoming entrained in the matrix flow and being subducted to depth and exiting the ROI.

Implications
The primary method of investigation is physical tectonic-fluid dynamic analogue laboratory modeling of the asthenospheric upper mantle. The primary variables explored are the style and geometry of subduction, the density contrast between the working and injected fluids, and the volumetric flux emanating from the source. This volumetric buoyancy flux can be correlated to a thickness of the unstable layer entering the subduction zone (see Figure 8 in Marsh (1979)). Analogue modeling has a long and rich history and has been used extensively in both geology and fluid dynamics [Hall, 1815;Hubbert, 1937;Riehl and Fultz, 1958]. See [Schellart and Strak, 2016]  The solid diapir experiments provided the ability to verify our fluid results. We do not want our external injection pressure to be too great and actually forcefully pump the diapir fluid vertically. This over-pressurization of the injected fluid is thought to produce more vertical paths because 1) it would maintain a connection with the source longer, allowing for diapirs to grow to be too large, and 2) we would be adding an upward vertical component to the diapir in addition to its own buoyancy.
Because the fixed-buoyancy solid and flux-buoyancy fluid diapir pathways look similar we assert that we are not forcing fluid into the system in an unphysical manner and looking at different physical processes. The desire is for our fluid injection to be as passive as possible and have the injection process itself minimally disturb the flow.
Essentially the goal is to push the fluid just hard enough to flow out from the source and immediately become entrained in the flow like dye being released and drawn passively into a stream. We are modeling solid state creep deformation processes of the mantle using a viscous fluid. Approximating the asthenosphere as a viscous fluid is acceptable over geologic time scales. It is then also appropriate to model the less viscous, buoyant diapirs with a fluid as well to better encapsulate the interactions we are seeking to explore. As discussed in our findings near-vertical ascent was not observed even with excessively buoyant material for the entire suite of subduction geometries, convergence rates, and buoyancy fluxes tested. Again, the density contrasts explored here are somewhat larger than expected for mantle wedge diapirs, which would be something on the order of 25%. However, if vertical source-vent pathways cannot be produced with the buoyancy that arises from our experimental density contrasts, it is unlikely that a smaller density contrast would result in diapirs experiencing less deflection. Therefore, we explored a range of possible density contrasts required to create near vertical ascent through the mantle wedge.
One of the buoyant fluids also contained water, which can have a strong effect on local matrix viscosity, and consequently diapir behavior in the analogue setting.
Water is thought to play a similar role in the natural mantle wedge prototype, lowering viscosity at the slab-wedge interface and within the actively convecting wedge facilitating rapid transport [van Keken, 2003;Grove et al., 2006;Wada and Wang, 2009;Jadamec, 2016]. Glucose rheology is heavily dependent on water content; there is roughly a 30% reduction in viscosity for each 1 wt% water added at 20 °C [Schellart, 2011]. The water in the wet ethanol diapirs can weaken the ambient matrix viscosity locally by up to two orders of magnitude. These localized weak zones are small features that become entrained in and advected by wedge flow, constantly redistributing heterogeneous features in space and time. Large scale chemical diffusion of this added water is negligible on experimental timescales, and the large volume of the tank relative to the ROI (ROI is less than 3.6% of the total tank volume) compensates for this effect. Hydrated ethanol diapirs are able to fully coalesce and reorganize into single bodies. This is thought to be attributable to water and its effect on bubble and droplet surface tension and viscosity at the material interface.
Natural diapirs themselves may melt by any of the three primary mechanisms: 1) external heating, 2) pressure release, 3) H20 flux altering the solidus. Sea-floor sediments and altered oceanic crust descend into the subduction zone at the surface of the cold lithospheric slab. The slab and accompanying material is cold relative to the surrounding mantle and as such diapirs are heated externally by the wedge throughout their ascent. The subducted diapir protomaterial may also be hydrated, resulting in a lower melting point than the surrounding wedge material. As the buoyant diapir material ascends, decompression and pressure release melting occur. All three of these are viable mechanisms supporting melt generation of and by diapirs as they traverse the mantle wedge, see Appendix 2. Water plays a key role in melt generation processes in subduction zones at and above the slab-wedge interface [Davies and Stevenson, 1992;Schmidt and Poli, 1998;van Keken, 2003;Grove et al., 2006;Hirschmann, 2006;Kelley et al., 2010]. Dehydration of hydrous phases at depth releases water that locally depresses the solidus initiating partial melting of the peridotitic mantle wedge.
If diapirs carry water with them into the wedge they will leave behind a trail of wet heterogeneity with enhanced melt potential. These trails are also weak allowing for subsequent diapirs to follow in these channels and ascend rapidly [Hall and Kincaid, 2001].
There is also some variance in constraining the rheologic flow law for olivine and bulk mantle compositions [Karato, 2010;King, 2016]. There are certainly many unknown or poorly constrained parameters that are of great importance to the overall dynamics of subduction zones and the mantle wedge. However the past and present plate rates that are supplying material into modern subduction zones are well known, and we leverage this as a starting point to investigate how the geologic protomaterial for rocks we observe at the surface interacts with the large-scale wedge flow.
Resolving all of the complexities of the subduction system is computationally expensive and complex. There are many dynamic feedbacks to the subduction system that are poorly constrained due to lack of directly measurable evidence. Until experimental and remote sensing efforts improve our ability to probe subduction zones and better constrain the system we must rely on theoretical, numerical, and analog models to gain insight. Altering some of the most basic assumptions commonly made about the subduction system can have drastic effects on model results and ultimate interpretation and understanding. The composition of the mantle wedge is not known with great certainty and is often assumed to be homogeneous with a bulk peridotite mineralogy dominated by olivine. Actual mantle xenoliths and ophiolite sections provide evidence that the sublithospheric mantle and wedge is not well mixed and can be highly heterogeneous in the compositional, thermal, and stress fields in both space and time [Behn and Kelemen, 2006;Behn et al., 2007;Le Voyer et al., 2017]. Global Rayleigh numbers for mantle convection are 5e5 and 5e7 for layered and whole mantle convection respectively [Turcotte and Schubert, 2014] [Gerya and Yuen, 2003;Hasenclever et al., 2011]. Three space dimensions and time evolution are necessary to evaluate this complex system, especially if we wish to evaluate trench-parallel transport along strike.

Conclusions
A suite of analogue models was used to examine the interactions between        The saddle point with no apparent motion relative to the trench is very clear at the right of the image (star shape). This point migrates with the trench and demonstrates why frame of reference is so important when considering rollback style subduction. Material that appears to be translating to the left is actually still moving to the right slower than the trench migration rate.      Tables 1,2) . a) The initial placement of rigid sphere diapirs is shown, where 4 diapirs are distributed in a trench-parallel orientation along each slab segment. As in prior cases, diapirs are initially placed in the VBL above the subducting plate. Colored circles are used to construct distinct diapir pathlines for those initiating above the steep belt (b) and shallow belt (c). (b) Patterns clearly show that flow/diapir transport is from the deep wedge above the steep slab segment to the sub-arc region above the inner edge of the shallow slab segment. c) Without rollback, diapirs sourced from the shallower slab also return to the sub-arc of this plate through curved ascent paths (1. slightly up and away from the trench, 2. near-vertical through midwedge, 3. strongly deflected sub-horizontal and toward the trench in shallow wedge). Grey shaded circles represent starting points along diapir particle paths.   Figure 15, but for case E_14, with rollback subduction and diapirs of (dry) ethanol. a) Frame shows a diapir that has risen from VBL into the mid-wedge, with a trailing conduit that has severed from source region and two smaller unsuccessful diapirs being entrained out of the ROI. b) As the diapir rises through the mid-wedge, the tail is drawn downward with background flow and is severed from the source.  Figure 18. Cross-sections with diapir paths for case E_7 as in Figure 14. (ethanol-dry diapir, downdip-only slab, 45° dip, 1 diapir source) highlighting a second mode of evolution/ascent. Paths are shown in blue for three distinct Mode II diapirs. These upwellings have received buoyant fluid from both up-slab and down-slab conduits while in the VBL, achieving a larger volume and rise velocity. Surfacing locations arẽ 150 km closer to the trench than Mode I (green paths shown for reference). A remnant down-slab conduit is seen as a bright tail in (b), trailing the second high volume (Mode II) diapir to surface in this E7 experiment. Figure 19. a) Cross-section through wedge for case E7 as in figures 17 and 18 highlighting the third mode of diapir interaction/ascent from slab to surface (shown as red circle paths). a) A repeatable pattern in these cases is that small diapirs separate more slowly from the VBL, and are dragged deeper into the system. Single diapirs stall here, where upward buoyancy balances downward drag. b) Accumulated smaller diapirs rise together at a greater ascent rate, along nearly vertical paths. c) Diapir clumps surface at the greatest distance from the trench, ~75 km further from the trench than reference (green Mode I) non-interacting diapirs. d) An interesting (repeatable) result is that inner diapirs in the clump can be delayed in surfacing (yellow path) by being held beneath their partners. . Frames cropped to focus on the southern slab segment. a) In 3-D experiment, even simple downdip-only subduction generates wedge flow (above the VBL) with a trench-parallel component that carries material from the slab edge towards the wedge centerline (see Figure 9). b) Diapirs moving in the trench-parallel direction interact with new buoyancy sources. c) In the 5 second (167 kyr) interval between frames the deeper diapir from source 4 merges upward into the shallower diapir from source 3. Volume and rise rate increase markedly and ascent becomes near vertical.

Figure 21.
Map-view images from experiment V_7 (downdip-only subduction), using vodka as the diapir fluid, that allows for merging of distinct flow structures. a) Trenchparallel drift of diapirs from different source locations produces a vertically aligned column of diapirs originating from three different source locations. b) In the 5 seconds (167 kyr) between images, the deepest diapir from source 4 has merged upwards into the mid-level diapir from source 3. c) In the 5 second (167 kyr) interval between subsequent frames the newly enlarged combined mid-level diapir accelerates and merges upward into the diapir from source 2 making one high-volume nearsurface diapir containing material from sources 4, 3, and 2. Figure 22. Side-view, trench-normal cross-sections through the wedge for experiment V_7 (downdip subduction, 4 vodka sources). a) Three distinct diapirs have moved into the mid-wedge, where they drift with a trench-parallel component of motion (into the page). b) One time-step later a deeper diapir from a different source is moving upward and into the page through the tail of the overriding diapir. c) Five seconds further in time the two diapirs have fully merged.

Figure 23.
Side-view, trench-normal cross-sections through the wedge for experiment V_7 (downdip subduction, 4 vodka sources). a) An array of diapirs are evolving at different trench-parallel locations (into the page). b) Deep diapir fluid encounters a conduit structure and begins rapid ascent. c) Five seconds ( 167 kyr) further in experiment time, the deeper is merging upward into the shallower diapir originating from one source closer to the centerline. c) After another five second (167 kyr) timestep, the newly merged diapir is beginning to merge upwards into the stalled diapir above. The deeper diapir fluid from (b) transits ~200km in 0.4 Ma, 12 seconds of lab time, yielding an average rise rate of 50 km / 100 kyr.  Map-view images from experiment V_8 (rollback subduction, 35° dip, 4 vodka diapir fluid sources, free surface boundary). The interaction scenarios are a complex balance of downward entrainment in the VBL and advection with the slab, vertical rise due to buoyancy flux and trench-parallel drift due to three-dimensionality of mid-to shallow-wedge return flow towards the sub-arc wedge apex and slab centerline. a) The middle of the three established diapirs is still connected to source 1, growing and taking in fluid from below. b-c) Three distinct diapirs are fully detached from their respective sources and are advected towards the wedge apex without further interacting with deeper diapir fluid. Figure 26. Side-view, trench-normal cross-sections through the wedge for experiment V_8 (rollback subduction, 4 vodka sources, free surface). Example of shallower diapir moving with strong trench-parallel motion and interacting with multiple, deeper diapirs arrayed along-strike. a) Mid-level diapir (outlined) is moving into page. b) (15 seconds elapsed from frame a) Deeper, growing diapir in VBL moves rapidly upward into overlying diapir. c) (15 seconds elapsed from frame b) The mid-level feature has moved in a trench-parallel direction over the next deep growing diapir (outlined). d) (5 seconds elapsed from frame c) Deeper buoyant fluid has begun rapid upward assimilation into overlying feature. e) Large volume diapir is impacting underside of plate after accumulating deep fluid from multiple sources in ~200 km of along-trench motion. Figure 27. Side-view, trench-normal cross-sections through the wedge for experiment V_10 (rollback subduction, 4 vodka sources, no-slip overriding, back-arc plate). a) An array of diapirs are evolving at different trench-parallel locations. The presence of a no-slip condition under the back-arc plate causes return flow pathways to deepen in the wedge. b-c) Deep diapir fluid encounters a conduit structure and begins ascent, assimilating into overlying diapir. d) The large volume diapir develops a larger ascent velocity. An enhanced trailing conduit provides a drifting pathway for much deeper fluid to ascend to the surface rapidly (~200 km in ~0.3 Ma (5-10 seconds)).

Figure 28.
Multiple still frames from experiment D_14, left panels are top view and right panels are side view. The diapir highlighted in magenta moved very slowly during this part of the experiment translating and rising only ~100 km over 15.3 Myr shown here, while its neighbors highlighted in blue and orange were advected much farther. This is because of the interaction of the magenta diapir and the pocket of material being extruded from between the opening slab tear (white dashed region).

Figure 29.
Map-view plot of along-strike transport paths and diapir positions relative to the trench. Markers are colored and sized by time. Spacing between markers is indicative of translational velocity (farther = faster). Black stars mark initial positions of rigid spheres at the slab-wedge interface, cyan stars mark position halfway along path, and red stars indicate final surfacing position. Note greater velocities at depth while diapirs are proximal to the kinematic forcing of the slab, and increased alongstrike displacements for diapirs originating farther from the slab centerline. 4 diapirs surfaced along the centerline but took very different paths to get there.    highlighting the regularity in air diapir timing and spacing. These diapirs form a nearly continuous chain from source to surface occupying the entire mantle wedge, giving a sense for the displacement these diapirs experience in each of the three wedge flow regimes. Diapirs in this experiment surfaced nearly vertically above the source, but were still deflected ~70 km away from the trench before entering the wedge apex return flow.

APPENDIX 2 Diapir Thermal Model
A key aspect to mass transfer through subduction zones is the thermal evolution of the material as it moves from the slab surface to the Earth's surface. In We chose a subset of observed laboratory diapir paths to force the boundaries of a 2D diapir heat model. The results from two end-member modeling runs are presented and compared, a short fast diapir path and a long slow diapir path. Final states for the two diapir end members are shown in Figure A2.5 with their accompanying temperature-time curves in Figure A2.3. All modeling runs are initialized with matrix and diapir temperatures equal to the slab surface temperature.
When internal flow is included heat is mixed around the inside of the diapir and the central temperature changes faster than for conduction only cases ( Figure A2.6).
Regardless of whether internal flow is included both the short-fast and long-slow diapirs achieve internal temperatures in excess of 1200 degrees, hotter than the wet peridotite solidus and melting is expected to initiate before arriving at the LAB.
A 2D advection-diffusion model of a spherical diapir with axial symmetry is forced at the boundaries by the lab-derived vertical velocity-temperature-time pathways to inspect the internal thermal evolution of an initially cold chemically buoyant diapir. This modeling is evaluated using 2D thermal models because 3D wedge temperature fields for experimental subduction geometries and convergence rates are not known to this author. This numerical model solves a dimensional temperature equation (eq. A1) for the thermal field using an explicit forward Euler finite difference formulation.  [Hadamard, 1911] and [Rybczyński, 1911]. The model domain is 20 [km] x 20 [km] with 165 x 165 standard A-grid nodes [Arakawa and Lamb, 1977] that solve for temperature after operator splitting is employed to first advect and then diffuse the temperature field at each time step. The internal flow is driven solely by the coupling of the internal and external fluids at the material interface and scales linearly with the external free stream velocity, so is scaled by the known vertical velocity from the lab trajectories. Numerical experiments were computed in two ways: 1) with internal flow and advection of heat included and 2) with conductive heat transfer only and no internal circulation.
This modeling effort helped gain further understanding at first order of the internal temperature structure and thermal evolution of a diapir and its potential for melting along the observed trajectories. Diapirs traversing the wedge gain enough heat from the surrounding mantle to raise their internal temperature above the wet peridotite solidus for both solid conductive heating and internal flow calculations. We employed a novel approach to modeling the thermal evolution of compositionally buoyant diapirs numerically. This effort shows that even large or fast diapirs absorb an appreciable amount of heat as they transit the wedge and internal temperatures surpass the wet peridotite solidus before arriving at the LAB.       This image highlights the nature of the internal circulation to transport heat from the diapir/mantle interface into the diapir.