The mechanical behavior of the crust and upper mantle materials is temperature controlled; therefore, we first examine the geometry effects on the thermal structure.
Two competing factors affect the thermal evolution of subduction zones. The subducting slab cools its environment, but the forced return flow brings in hot material. Kincaid and Sacks (1997) conducted numerical convection experiments in a viscous, Boussinesq, Newtonian fluid within a two-dimensional Cartesian geometry. Fluid motion was due both to thermal convection and slab penetration. Slab geometry was not specified; it was found to be sensitive to slab buoyancy and slab velocity.
For parameters appropriate to northeast Japan subduction, a rather similar linear slab was found. The time-dependent evolution was examined for a wide range of parameter space. We rely on some conclusions that were found to be quite robust. For a relatively shallow slab dip angle such as is appropriate to present-day northeast Japan, cooling by the slab dominates. For steeply dipping slabs such as the Marianas, the hot return flow controls the temperature.
Importantly for our study, the shallow dip case leads to a runaway process, which dominates the mechanical interaction of the subducting plate and the overlying structure. Because at this geometry the cooling of mantle material near the apex of the mantle wedge both to the overlying crust and to the slab overwhelms the heating by the return flow, the material viscosity increases and it tends to be more sluggish. This gives it more time to cool, and so its viscosity continues to increase until it becomes almost static. Eventually, this frozen "nose" extends from the apex and nearly reaches the volcanic front. The time required to reach this condition is estimated to be on the order of 3 m.y. Figure F11 shows the nose growth as a function of time. The effect of this now rigid wedge is to increase the coupled length between the slab and the overlying plate.
Honda (1985) also constructed a model of thermal structure of the northeast Japan area to fit the heat flow observation using available geophysical, geochemical, and geological constraints. His modeling required the existence of a rigid (nonmoving) mantle nose region with low temperature bounded by the subducting slab extending about halfway to the VF from the trench. This portion corresponds to the high-velocity mantle area (VP = 8.0 km/s) (Fig. F3).
Mechanical coupling affects mountain growth and continental deformation. Huang et al. (1998) examined the deformation of northeast Japan resulting from long-term (5 m.y.) interaction with the Pacific slab. A finite element model, TECTON, from Melosh and Raefsky (1980) and Wallace and Melosh (1994) was modified to include dynamically applied fracture criteria and erosion and deposition. Because the strike of folds and faults in northeast Japan are approximately normal to the direction of subduction, a two-dimensional plane-strain analysis is adequate. The strength profile of the lithosphere was based on experimentally determined rheology laws (Wilks and Carter, 1990) and temperature derived from heat flow data (Furukawa and Uyeda, 1989). The criteria for a successful model are that the topography, gravity residual anomaly, and seismicity must match the observations.
The output of the model is the coupling force between the subducting slab and the overlying continent and its depth distribution. Additionally, the erosion and sedimentation rates are determined. The huge gravity residual anomaly of +75 mGal over the Kitakami range is supported by the warp generated by the strong coupling of the Pacific plate down to ~50 km. This is consistent with the observation of thrust earthquakes on the interface down to this depth (e.g., Hasegawa et al., 1994). The mantle wedge must be rigid enough to not only store strain released by earthquakes but also to provide long-term support for the warps.
The growth rate of the near coastal range, the Kitakami, as well as the Backbone range is proportional to the compressional force applied by the subducting Pacific plate. Ohmori (1978) showed that the erosion rate is proportional to the square of the elevation, and his estimated value of 0.3 mm/yr for a height of 1000 m was used in the model runs (Huang et al., 1998). Whereas erosion is readily calculated, deposition of the erosion products is more complicated. Downwarps on land get filled, but we assume much of the sedimentation can be represented by offshore depositions.
The results from the two models described above can be used to relate the proposed changes in the subduction geometry to observables such as sedimentation rate and volcanic flux determined from the ODP and DSDP holes in the forearc (Figs. F6, F7).
The proposed model geometry showing plate-coupled zone length spanning 20 m.y. is shown in Figure F12. At 20 Ma the slab geometry is very similar to that of present-day Marianas. The distance of the postulated volcano identified in Hole 439 to the trench is similar to that in the Marianas arc. Backarc spreading occurs in this arc, and by analogy with Japan at ~20 Ma, tension in Japan and opening of the Japan Sea would be expected.
By 14 Ma the volcanic front reached the Japanese island, shown in Figure F9. As well as the dip change, there also had to be tectonic erosion of the overlying crust starting at least at 20 Ma. This erosion was continuous, as the graben formed by the bending of the Pacific plate could scrape material off the bottom of the overlying plate. This caused the slab and, therefore, the volcanic front to advance by ~30 km/10 m.y. (von Huene et al., 1994). In addition, at this stage the slab may not have been at its equilibrium shallow dip yet. This also caused westward migration (Fig. F10). The shallow dip angle encouraged cooling and eventual freezing out of the mantle wedge. By 6 Ma this stiff zone reached its maximum (equilibrium) extent and the mechanical interaction between the subducted slab and the overlying material was also at a maximum.
We have so far only discussed the thermal effects of the subducted cool slab on the temperature-dependent viscosity of the overlying mantle. In addition, water is released from the slab. Ultimately, some of this water flux induces melting in the mantle and causes volcanism. Recent measurements of boron isotopes suggest that the slab releases 11B-enriched fluids from the shallowest levels to depths below 200 km (Benton et al., 2001). If this water is released into material with a temperature below its wet solidus, hydrous phases will form. Eventually, the static mantle wedge will become hydrated. This causes lower seismic velocities (Fig. F3) and intermediate resistivity. Later, the vicinity of the slab interface will become saturated and further water release will give rise to high pore pressure and, finally, free water. Since this is a cumulative process downslab, the excess will first appear at the greatest depth (i.e., the deepest part of the static wedge).
The effect of the high pore pressure is to reduce the mechanical coupling (shear) between the slab and the overlying static mantle. This reduced coupling length has a strong effect on the eastern mountain building. As shown in Huang et al. (1998), the Kitakami range requires the force couple provided by the deeper locked zone. Reduction of the coupling length affects the mountain building and therefore the erosion, which is the only observation, in two ways; both the coupling torque and the coupling length are reduced.
In order to model the effect of the coupled length on the resultant force, we follow the reasoning of Wang and Suyehiro (1999). Making the simplest assumption that friction was constant, it follows that the force on the overlying plate would be simply proportional to the coupled length. The coupling lengths are shown schematically in Figure F12. The slab contact length in the crust depends on erosion, as is shown in Figure F10.
Using water depth as a determinant, erosion has been modeled by Lallemand et al. (1992) and von Huene and Lallemand (1990). Since the water depth is affected not only by the erosion but also by the sediment load and the slab coupling, the calculation is complex. For the purpose of this model, we assume a simple progression because we only require a plausible rate. We applied the modeling techniques described above to the slab geometry changes postulated from the volcanic arc position shown in Figure F10. The primary data that we had to match were the sedimentation rate and volcanic flux (ash) data determined from the cores taken during ODP Leg 186. Figure F6 shows the parameters of the calculation. The results are shown in Figure F13. Although it is not possible to determine the deposition paths well enough to quantify the sedimentation offshore, we can calculate the profile of the temporal changes in sedimentation rate.
The effect of water saturation leading to high pore pressures will not only reduce the coupling length and limit the seismic thrust zone depth (to ~55 km) but will also affect the volcanic flux. The wedge had cooled sufficiently to become essentially static and to couple strongly by ~6 Ma, though the coupling may have started increasing from ~10 Ma (Figs. F10, F12). The overlying mantle became increasingly hydrated and less absorbent until ~4 Ma, when high pore pressure and free water existed below the seismic thrust zone. The free water traveled with the slab until it reached hot, mobile mantle with a temperature at least above its wet solidus; there it will flux melt. This water added to that which comes from dehydration of slab minerals such as amphibole, so the volcanic flux increased. Note that the model prediction is that the increase in volcanic flux should occur when the sedimentation rate decreases, as is observed.