ZHANG Hao, GONG Wei, 2), XING Junhui, 2), XU Chong, and LI Chaoyang

The Subduction Structure Beneath the New Britain Island Arc and the Adjacent Region from Double-Difference Tomography

ZHANG Hao1), GONG Wei1), 2), XING Junhui1), 2), XU Chong1), and LI Chaoyang3), *

1),,,,266100,2),,266237,3),,266590,

We applied double-difference tomography to relocate seismic events and determine the lithospheric velocity structure be- neath the New Britain Island arc and the South Bismarck Sea Basin, based on the local P wave arrival time dataset collected by the International Seismological Centre. Results of the seismic relocation and velocity inversion show that the subduction of Solomon Sea Plate along the New Britain Trench is spatially different above 150km, and the subduction angle of the slab on the west side is higher than that on the east side. The relocated earthquakes also show that there are double seismic zones at the depths of about 30–90km beneath the New Britain Island Arc. The velocity structure shows that the dehydration of the subducting slab caused the low-velocity ano- malies in mantle wedge above the slab, which are associated with the magmatic activities around the New Guinea-New Britain Island arc. Moreover, it shows that there is another low-velocity anomaly zone beneath the Bismarck mid-oceanic ridge with spatial variation. Beneath the west of the Bismarck mid-oceanic ridge, the low-velocity anomaly is weakly connected to the subducted Solomon Sea slab. Conversely, the low-velocity anomaly beneath the Manus Sea Basin is highly intertwined to the subducting slab and its mantle wedge, indicating that the subduction of the Solomon Sea Plate might be a key deep dynamic factor that drives the spreading of the Manus Sea Basin and the separation of the Bismarck Plate.

New Britain Trench; Bismarck Sea Basin; Manus Sea Basin; double-difference tomography

1 Introduction

The New Britain Island Arc is an essential part of the New Guinea-Solomon arc system in the Southwest Pacific, which is located in the convergent boundary between the Pacific Plate and the Indo-Australian Plate. Multiple tectonic activities caused the creation of complicated subduction systems, including a large number of ‘trench-arc-basin’systems here (Fig.1). The complex subduction patterns, fre- quent seismic events, and active magmatic volcanic activities make this region become one of the most active areas with intense crustal movements and complex geological processes among all the subduction zones (Gong., 2019; Feng., 2022). Meanwhile, with those characteristics, it is also an ideal region for investigating the subduction structures and dynamics.

The Bismarck Basin, located on the north side of the New Britain Island Arc, is a typical back-arc spreading basin throughout multiple complex tectonic periods. About 10Myr ago, the subduction of the Pacific Plate at the Manus Trench was terminated because of a collision with the Ontong-Java Plateau. The collision reversed the subduction and created a new subduction zone at the present location of the New Britain Trench, which caused the northward subduction of the Solomon Sea Plate (Lee and Ruellan, 2006; Franz and Romer, 2010). About 6Myr ago, the Bismarck Sea began to expand as a back-arc basin, andthen, it was transformed by the later arc-continent collision (Holm., 2016). At 3.5–4Myr ago, the New Britain Island arc-Huon Peninsula collided with the New Guinea continent, making it rotate counterclockwise, and then, the Bismarck Basin began to split, forming the Manus Basin (Tregoning., 1998; Weiler and Coe, 2000; Holm., 2016). At the same time, the rotation of the collision zone between the New Britain-Huon Peninsula and the New Guinea continent formed three large left-lateral transform faults (Willaumez, Djaul, and Weitin transform faults). The active center of the basin is now concentrated on the Manus Basin, which has been spreading asymmetrically for 3.5Myr (., Taylor, 1979; Both., 1986; Binns and Scott, 1993; Zhao., 2014; Ma., 2017). The mid- oceanic ridge and transform faults in the Bismarck Sea, the Ramu-Markham fault, the Manus Trench, and the NewBritain Trench divide the Bismarck Plate into multiple micro-blocks, which include the North Bismarck block, the South Bismarck block, and the Manus micro-block.

Massive earthquakes have occurred in the lithosphere of this area, triggered by volcano activities, transform faults, mid-oceanic ridge, and subduction. These intense seismic activities have provided the possibilities to study the spatial and structural characteristics of the micro-blocks and subduction zones. Based on the analysis of crustal velocities derived from the deep seismic reflection data, the Moho surface was found to deepen gradually from the Bismarck Sea toward the New Britain Island (Finlayson., 1972).Taylor (1979) demarcated the boundaries between the sou- thern and northern blocks of the Bismarck Basin according to the seafloor topography, focal mechanism, and seafloor magnetic anomaly data and presented an overview of the basin structural characteristics. The data of natural epicenter and focal mechanism demonstrate that the subduction angle of the Solomon Basin along the New Bri- tain Trench reaches about 70˚, whereas the dip angle of the slab is less than 40˚ at the Trobriand Trench (Abers and Roecker, 1991; Pegler., 1995). Moreover, the distribution of epicenters was utilized to determine the subduction pattern of the slab along the New Britain Trench. Also,on the basis of this data, plate reconstruction has been conducted to restore the tectonic movements in this area since 8Myr ago (Holm and Richards, 2013; Holm., 2016).

Fig.1 Structural map showing the tectonic setting and geological tectonic units. The top-left inset is a geographic map of the West Pacific and surrounding regions, in which the blue rectangle shows the location of the study area. The red lines denote the vertical sections with the tomography results in Fig.10. DF, Djaul fault; HFB, Huon-Finisterre block; MNB, Manus block; MT, Manus Trench; NBI, New Britain Island; NBS, North Bismarck Sea; NBT, New Britain Trench; NGI, New Guinea Island; NGT, New Guinea Trench; RMF, Ramu-Markham fault; SBS, South Bismarck Sea; SS, Solomon Sea; TT, Trobriand Trench; WIF, Willaumez fault; WTF, Weitin fault.

From the above, the accurate seismic location and lithospheric velocity structure have still been scarce in this area. For this reason, we collected P wave travel-time dataset in this region and applied double-difference tomography (Zhang and Thurber, 2003) to obtain the more accurate distribution of the earthquake epicenters and the more reliable 3D lithospheric velocity structure. Combined our tomography results with previous studies in this region, we seek to cha- racterize the subduction geometry and the mantle structures in the lithosphere, identify the melting mantle wedge induced by the subduction, and provide seismic tomographyevidences for understanding how the New Britain subduction zone affects the Bismarck Basin.

2 Data and Method

2.1 Data

The data used in this study is derived from the current catalog dataset provided by the International Seismological Centre (ISC). The range of study area is between 142˚–156˚E and 2˚–10˚S. There are a total of 85 seismic stations deployed in this area, which are published on the ISC website. We retrieved the absolute P wave arrival time data of the local earthquakes that occurred from January, 1960 to January, 2017, in the area with slightly different range between 142˚–156˚E and 2˚–10˚S. To improve data quality and ensure a reliable seismic inversion result, we filtered the original data by removing abnormal travel timeinformation using the time-distance curve method. The filtered data was used as the input for double-difference tomography, that is, the travel-time data between the two green lines in Fig.2. In Fig.2, most travel time data from ISC converge near the fitted time-distance curve, which also indicates that the time information in the ISC databasehas been selected in advance and corrected by relevant staffs and is therefore relatively accurate and reliable. Double-difference tomography requires event pairs to performseismic relocation and velocity inversion. Therefore, we set the following parameters in our study: the maximum dis- tance between two stations of event pairs is 800km, and the minimum and maximum distances between the source locations of event pairs are 0.1 and 50km. Each event can be combined with up to 50 events to form event pairs. The minimum and maximum numbers of phases required for each event pair are 3 and 50, respectively. After setting theseparameters, we obtained 6659 suitable seismic events, 32352 absolute arrival data, and 383086 relative arrival data for later inversion. The distribution of the stations, events, and ray paths is shown in Fig.3, and the ray paths cover most of the inversion region, allowing to extract relatively reliable inversion results.

Fig.2 Diagram for fitting travel time with the epicentral dis- tance.

Fig.3 Distribution of ray paths in the study area. The red dots indicate the source locations, the black triangles indicate the station locations, the blue plus symbols indicate the inversion grid nodes, and the gray straight lines indicate the ray paths for tomography.

2.2 Double-Difference Tomography Method

Zhang and Thurber (2003, 2006) proposed the double- difference tomography on the basis of the double-diffe- rence location method (Waldhauser and Ellsworth, 2000), which is an improved version of travel time tomography. The theory of double-difference tomography is as follows:

According to the ray theory, the earthquake travel time can be expressed as Eq. (1), whereis the observed ar- rival time of the body wave from the sourceto the station,is the time of earthquake occurrence,is the slowness vector, and ds is the integral path element.

For local earthquake tomography, the residualbetweenthe theoretical arrival time and the actual observation time from the seismic eventto stationcan be expressed as:

In Eq. (2), the location of the seismic source (∆1, ∆2, ∆3), the time of occurrence, the slowness field, and the ray path are all unknown. If another eventis recorded by the same station, an event pair can be formed. The differences of travel time residuals of these two events are double differences:

Double-difference tomography uses both the absolute and the relative arrival time data for the joint inversion of the three-dimensional velocity structure and the source lo- cation. Since the location of the seismic source affects the inversion accuracy of the velocity structure in the local seismic tomography, the first few iterations of the actual inversion will set a larger weight value for the absolute time data to establish the rough three-dimensional velocity structure of the entire study area. Then, the relative arrival time data will be given a larger weight value in the next iterations to constrain the velocity structure near the sourcearea. The velocity structure inversion and the seismic relocation are calculated simultaneously, so the velocity structure can be continuously inverted and adjusted to obtain a higher accuracy of the source location results compared with the double-difference location method. Moreover, unlike the traditional travel time tomography method, relatively accurate seismic location results and the use of relativetime data allow to obtain a more reliable velocity structure (Zhang and Thurber, 2003).

The regional double-difference tomography (tomoDD) used in this study applied the ‘node method’ to parame- terize the model. The ‘pseudo-curve method’ is used to cal- culate the travel time and ray path. The inversion problem was solved by the damped least squares algorithm (LSQR, Paige and Saunders, 1982). The P wave double-difference travel time residual data and absolute travel time residual data in the study region are used to obtain more accurate results of the source location and three-dimensional velo- city structure.

2.3 Model and Parameters

The range of the inversion area is 143˚–155˚E and 2˚–8˚S. The inversion area includes the New Britain Island, the South Bismarck Basin and its adjacent region. The grid interval is 0.5×0.5. The grid nodes are shown as the blue plus symbols in Fig.3. Considering the focal depth of the local earthquakes published by ISC, the node depth in the vertical direction is set as 10, 30, 50, 70, 90, 120, and 150km. Since this area lacks a relatively accurate velocity model, the initial velocity model used in the inversion is the adjusted velocity model of IASP91. The initial velocity at each depth node is listed in Table 1.

Double-difference tomography applied LSQR to iteratively solve the inversion problem. During the inversion process, the damping and the smoothing parameters control the convergence speed of the inversion and the smoo- thness of the constraint model, which have a great impact on the reliability of the inversion results. Thus, a suitable damping parameter and a smoothing parameter must be selected before the inversion. The best parameter combination should make the model variance change slightly and the variance of the inversion result decrease quickly. In this study, the L curve method is used to test the damping and smoothing parameters. By comparing the inversion re-sults based on different parameter values, the trade-off cur- ve between residual norm and model norm was obtained, and the corner point of L-shaped curves represents the va- lue of optimum balance between time residual minimization and velocity model smoothness. The L curves of two pa- rameters,., smoothing and damping, are shown in Fig.4 and the optimal values are selected to be 20 and 100.

Table 1 Initial model of P wave velocity

Fig.4 Optimum smoothing parameter (a) and damping parameter (b) selected by using the L curve method and marked by red circles.

2.4 Checkerboard Resolution Test

Before inversion, we applied the checkerboard test (Spak-man., 1993) to analyze the model resolution, which isa common method in tomography to check the quality andreliability of inversion results. The checkerboard test firstly generates a theoretical chessboard model by adding positive and negative velocity perturbation to the initial velo- city model. Then, the theoretical chessboard model and the real distribution of sources and stations are used to calculate the theoretical travel time as the test data. Finally, the model for analysis is inverted based on the test data. The model resolution is determined by comparing the inversion result to the theoretical chessboard model. If the checker- board model can be recovered in some areas, then we assume that the inversion results in this area are reliable. Forour case, the theoretical chessboard model is generated basedon the initial model with±5% perturbations and 0.5˚×0.5˚ grid intervals.

Fig.5 shows the results of the checkerboard test at each inversion depth. At the depth of 10km, the resolution in our study area is relatively low. Only in few parts, near the east of New Britain Island and the New Guinea Peninsula, the velocity perturbation was recovered because the seismic station distribution severely limited the shallow ray path. As the depth increases, the ray path gradually covers most of the inversion region. The inversion results at depths of 50,70, 90, and 120km have higher resolutions, suggesting that the model has a high resolution at the upper mantle.

3 Results

The velocity structure inversion and source relocation are simultaneously conducted in double-difference tomogra-phy, in which the relative travel time data of seismic eventpairs are used to improve the accuracy of the relative location between sources. The horizontal distribution patterns of earthquakes before and after relocation are shown in Fig.6, from which it can be seen that the changes of the ho- rizontal location are not obvious. Fig.7 shows the changes of focal depth before and after relocation. The number of earthquakes at depths of 0–20km is relatively reduced, whereas the number of earthquakes increases at depths of 20–40km. At other deeper depths, there are slight changes. Overall, the change of the epicenter location after relocation is extremely slight, which also reflects that the seismic information of the ISC database is relatively accurate and reliable.

Fig.5 Results of the checkboard resolution tests.

Fig.6 Spatial distribution of earthquakes before (a) and after (b) the relocation calculation.

Fig.7 Focal depth changes before (a) and after (b) the relocation calculation.

For the inversion results of velocity structure from tomo- graphy, Fig.8 shows the distribution of travel time residuals before and after the inversion. The range of time resi- duals drops dramatically from −5–5s to −0.1–0.1s after the inversion. Fig.9 shows the velocity horizontal slices ofboth the velocity structure and earthquake location at depths of 10, 30, 50, 70, 90, 120, and 150km.

For the shallow depths above 30km, due to the limit of resolution, we only analyze the areas where velocity pertur- bations are recovered or partially recovered in the checker- board test. As shown in Fig.9, there are obvious low-ve- locity anomalies beneath the New Britain Island arc, the New Guinea Peninsula and the Solomon Island arc. The Solomon Sea Basin and the Bismarck Sea Basin are cha- racterized by high-velocity anomaly zones. In the depths of 10–70km, there are belt-like low-velocity anomalies on the southern boundary of the South Bismarck Block, from the Ramu-Markham fault to the New Guinea Trench. The distribution of these belt-like low-velocity anomalies is consistent with the seismic zone along the boundary. As the depth increases, from the depth of 50km, the earthquake obviously converges to the vicinity of New Guinea-New Britain-Solomon Island arc from the west to the east; but the seismic zone along New Britain Trench gradually movedtoward the Bismarck Basin. Taking the seismic zone as the boundary of plate, it can be seen that the velocity of subducting slab on the southwest or south sides of the South Bismarck Basin is relatively high. This high-velocity ano- maly zone also has a good correspondence to the shape of New Britain Trench. At depths of 50km, an obvious low- velocity anomaly appears at the boundary of the North- South Bismarck Basin, close to the Manus Basin. Furthermore, a large low-velocity anomaly appears in the mid- west part of the basin at the depth of 90km. This anomaly can be tracked to the depth of 150km and reaches the ma- ximum size at the depth of 120km. Despite of the distribution variation at different depths, the velocity beneath southern Bismarck Basin is relatively low compared with the subducting slab on the boundary at the depth greater than 50km.

Fig.8 Distribution diagrams of the time residuals associated with the initial model (a) and ­final model (b).

Fig.9 Horizontal slices of the tomography results. The white lines denote the boundaries of plates.

4 Analysis and Discussion

Based on the relocated earthquakes and the velocity stru- cture, we extracted four vertical sections along the red lines in Fig.1 to discuss the velocity structure of the lithosphere beneath the South Bismarck Basin and the subduction pattern of the Solomon Sea Plate along the New British Trench.

4.1 Earthquake Distribution

As shown in Figs.6 and 9, the distribution of relocated earthquakes clearly indicates the boundaries between different plates in the study area. At the depths above 30km, the earthquakes near the mid-oceanic ridge of the Bismarck Ocean Basin are continuously distributed along the east- west direction, indicating the Bismarck Sea Basin diving into South and North Bismarck microplates. Previous stu- dies (Tregoning., 1998, 1999; Wallace., 2004; Holm., 2013, 2016) assumed that the different characteristics of movement between the South and North Bismarck blocks lead to the division of the whole ocean basin. Nowadays the North Bismarck microplate is almostcompletely coupled with the Pacific Plate (Tregonging., 1998; Wallace., 2004), as there is no obvious relative movement between the North Bismarck Basin and the Pa- cific Plate. However, the South Bismarck Plate, which is affected multiply by the relative movement of the Pacific Plate and the Australian Plate as well as the subduction of Solomon Sea slab along the New Britain Trench, presents a relatively clockwise rotating motion (Tregoning., 1998; Weiler and Coe, 2000; Holm., 2016). According to previous studies of the focal mechanism and stress field in this area (Johnson and Molnar, 1972; Mori, 1989; Zoback, 1992; Tregoning., 1998; Heidbach., 2010), there are a large number of strike-slip faults at the boun- dary between the South and North Bismarck microplate. The development of this nearly EW trending strike-slip fault system also suggests the relative movement between the southern and northern sea basin. In Fig.9, the seismiczone of the westernmost strike-slip fault is distributed alongthe EW direction, with the length extending to 200km. Earthquakes along this strike-slip fault still occur at the depth of 50km, indicating that the active depth is up to thebase of lithosphere. Furthermore, there are still sparse seis- mic activities at the depth of 50km in the WIF and the east side of the Manus Basin. As a relatively young basin, the Manus Basin is still spreading at a high rate (up to 137mmyr−1), making it one of the fastest spreading basins(Taylor., 1979; Ortega-Osorio., 2010). The earth- quakes at the depth of 50km also reflect its strong mag- matic activities within the lithosphere, corresponding with the low-velocity anomaly at the depths greater than 50km. The deep focal depth of 50km and mantle low-velocity ano-maly suggest that the young Manus Basin is at early sprea- ding stage (Taylor, 1979; Tregoning., 1999; Holm., 2016) with the inadequate cooling of the lithosphere near the spreading mid-ocean ridge.

At the New BritainTrench, the earthquakes in Fig.9 clear- ly shows the distribution of the subduction zone. Meanwhile, the high-velocity anomaly extends northward as the depth increases, indicating that the Solomon Sea Plate is subducting beneath the South Bismarck Basin. The earthquakes in this area usually occur near the top of subducting slab. However, previous studies (McGuire and Wiens, 1995) found that there are double seismic zones (DSZ) with different stress states (horizontal, along-strike compression in the upper zone, downdip tension on the lower zone). The DSZ is generally caused by the dehydration em- brittlement inside the slab, furtherly resulting in fractures and earthquakes (Yamasaki and Seno, 2003). As shown in Figs.9 and 10, the seismic zone widens at depth of 30–90 km, which probably indicates the existence of DSZ. The red dashed rectangles in profiles P1–P4 (Fig.10) also show that earthquakes occur not only near the top surface but also inside of the subducting slab. Zhang and Wei (2008) analyzed all the DSZs within the subduction zones around the Pacific and found that the DSZs within young subduc- tion zones usually converge at depth above 150km. As shown in Fig.9 in this paper, the seismic zone along the New Britain Trench was significantly narrowed at the depth of 120km, suggesting the convergence of DSZ.

In addition, at the junction between the New Britain and the Solomon Islands, the extending direction of the seis- mic zone changes from SW-NE in the west to NW-SE in the east, which represents two different subduction directions of the Solomon Sea Plate. As shown in Fig.9, the seismic zone at depths greater than 120km was gradually separated at this junction, leaving a seismic blank zone. This observation is also in agreement with the study of Holm and Richards (2013), which suggested that subducting slab may have been torn in this area. However, limited by the resolution of the seismic velocity in the deep zone, the tearing of slab is poorly understood and requires further detailed studies to provide more evidence.

4.2 Subduction Pattern

As the black dashed lines shown in Fig.10, the subduction patterns are significantly different in the lithosphere of different regions, as we interpreted according to both theepicenters distribution and velocity structure. Although theslab subducted at a generally high angle along the New Bri- tain Trench, the dipping angle seems to increases from east to west. In the P1 and P2 profiles, the subducting slab is almost vertical, while the dipping angles in the P3 and P4 profiles are 70˚. The higher dipping angle on the west side may suggest that the western subduction zone was modified by the arc-continent collision (Holm., 2016).

Particularly, as shown in the P1 profile, the number of earthquakes that occurred at depths of 30–120km is larger than that of other profiles because of the intense tectonic activity in this arc-continent collision area. Moreover, the earthquakes are distributed on both sides of the trench, indicating that the Solomon Sea Plate forms an ‘inverted U-shaped’ doubly dipping structural pattern (Pegler., 1995; Whitmore., 1997). The two low-velocity ano- malies on the both side of the New Britain Trench represent two different blocks. As shown in the P1 section, the New Guinea Peninsula on the left belongs to the Australian Plate, and the Huon-Finisterre block on the right is a part of the South Bismarck Plate. The shallow earthquakes beneath those two blocks were caused by the doubly sub- duction of Solomon Sea slab, while the deep earthquakes (extend to the depth of 120km) might be associated with arc-continent collision. At the convergent margin of triple plates, the subducting Solomon Sea slab was strongly reformed, and there might be more fractures inside the slab. Thus, the distribution of earthquakes in P1 section was extremely dense.

Fig.10 Vertical sections of the tomography results along the profiles shown in Fig.1. The color bar of velocity is shown at the bottom. The black dashed lines show the subduction pattern of slab as we interpreted. L1, L2, and L3 denote three low- velocity anomalies. The red arrows show the potential upwelling flows of magmatic materials. The red dashed rectangles denote the double seismic zones. Black dots mark the location of earthquakes. Red triangles represent the active volcanoes.

4.3 Low-Velocity Anomalies and the Mantle Wedge

As shown in Fig.10, at the shallow depth above 30km, all profiles show obvious low-velocity anomalies beneath the volcanic islands arc, suggesting the influence of magmatism. At the depth of 90–150km, there is a low-velo- city anomaly L1 in the Midwest of the South Bismarck Basin, which is on the top of the subducting slab. We interpret it as the mantle wedge, which is an import structural unit of the whole subduction system. It is conventionally believed that the wedge is the deep source of magma beneaththevolcanicarc(.,Tatsumi.,1983,1986;Wy- llie, 1988; Tatsumi and Eggins, 1995; Wyss., 2001; Stern, 2002; Hyndman and Peacock, 2003). The magmatic materials are generated in the mantle wedge due to infiltration of aqueous fluids and partial melts derived from subducting slab, and then immigrate upward into the crust, triggering the back-arc volcanic activities. This low-velocity anomaly in the upper mantle can be clearly tracked to the depth of 150km in profile P1 and P2. Another low-velo- city anomaly L3 is imaged in all profiles, which is also located on the top of slab but at shallow depth, indicating another mantle melting structure at the depths of 30–90 km. With a high convergence rate of 9–11cmyr−1at the New Britain Trench (Honza., 1989; Tregoning., 1998; Wallace., 2004; DeMets, 2010), the Solomon Sea slab carry large quantities of H2O to the deep. As the depth increases, the slab releases some water into the overlying mantle wedge, causing partial melting and upwelling flow of mantle materials, as shown by red arrows in Fig.10. Our profiles demonstrate that the upper mantle in the mantle wedge supply the magma for the volcanic activities near the island arc.

Beneath the spreading mid-ocean ridges in the BismarckBasin and the Manus Basin, there are both earthquakes and low-velocity anomalies, as shown in profiles P1–P4 from west to east of Fig.10. For the old west spreading ridge of Bismarck Basin in profile P1, the relation between the low- velocity anomalies under the ridge and subducting slab is not clear, but it seemed that both are associated with the deeper source in the north. Beneath the new spreading center on the east, as shown in profiles P2 to P4, there is a significant low-velocity anomaly L2 at the depths of 50–90km, and the connection between L2 and the subduction of the Solomon Sea Plate is gradually enhanced from the west to east. In the P3 and P4, there are plenty of earthquakes under the Manus Basin, suggesting the intense mag- matic activities in the lithosphere. It is already known that the Manus Basin is caused by the transition of the spreading ridge and the clockwise rotation of the South BismarckPlate (Tregoning., 1998; Weiler and Coe, 2000; Holm., 2016), which is a ‘passive type’ of back-arc basin (Mu., 2019). However, the low-velocity anomaly L2 and intense earthquakes within the lithosphere suggest thatthe spreading of Manus Basin may be also affected by upwelling materials from deep mantle wedge. As shown in eastmost profile P4, the low-velocity anomaly L2 might be generated from the mantle wedge above the subducting Solomon Sea slab, and epicenters beneath Manus Basin are mostly distributed on the side close to subduction, suggesting close relations between magmatic activities and the subducting slab. The hydrothermal studies in the eastern Manus Basin (Thal, 2014; Ma., 2018) proved that the hydrothermal fluid on the spreading ridge is more influenced by volcanism and subduction compared with normal mature ocean mid-ridge, resulting in a magmatic and subduction-type fluid. The geochemical study of mag- ma from back-arc volcanic zones in the Manus Basin (Sinton., 2003) also suggests that the metasomatism and enrichment of subduction-related components has occurred in the mantle. We thus believe that the subduction of the Solomon Sea Plate along the New Britain should representanother important deep dynamic factor driving the spreading of the Manus Basin. Therefore, the Manus Basin is not totally caused by ‘passive’ back-arc spreading, but the combination of ‘passive’ and ‘active’ spreading affected by both the relative movement of plates and the upwelling magmatic flow from the mantle wedge above subducted slab. As for the relationship between the above dynamic factors, further studies are required.

5 Conclusions

In this study, the double-difference tomography is applied to obtain the seismic location and three-dimensional lithospheric velocity structure of the South Bismarck Basin and its adjacent areas based on the P waves travel time data provided by the ISC database. After the inversion, the travel time residuals were significantly reduced, and the seismic location accuracy was improved. The relocation result shows that the shallow earthquakes clearly define the boundaries of the South and North Bismarck Plates. The earthquake distribution at the New Guinea-New Britain Trench also clearly indicates the spatial distribution and subducting direction of the subduction zone.

Through the tomography results in our study areas, it can be seen that the subducted Solomon Sea slab is spatially different in its subduction patterns: the subduction angle within the lithosphere of the Solomon Sea Plate decreases from west to east along the New Britain Trench; the subducted slab on the west side is almost vertical because of the modified arc-continent collision, and a double-dipping structural pattern is formed at the convergence zone of the New Guinea Peninsula block and the Huon-Finisterre block. According to the epicenter relocation, double seismic zones may occur at a depth range of 30–90km, which is associated with the dehydration of the shallow subducted slab. Additionally, the velocity struc- ture shows that there are low-velocity anomalies above the subducted slab, which corresponds to the partial melting mantle wedge structure, providing material source for mag- matic activities near the island arc. Also, low-velocity ano- malies are imaged beneath the spreading ridge of Manus microplate. Profiles show that this mantle low-velocity anomaly is characterized by spatial heterogeneity: the mantle anomaly on the west side is barely affected by the subducted Solomon Sea Slab, whereas the low-velocity ano- maly on the east side shows a strong connection with the slab, which indicates that the subduction of the Solomon Sea Plate on the east side became another deep dynamic source for the spreading of the sea basin. It can be seen from the inverted velocity structure and distribution of ear- thquakes that the subducted Solomon Sea Plate and its de- hydration played a key role in the spreading of the Manus microplate and the separation of the North and South Bismarck Plates, as well as regional volcanic activities. In the future studies, more geological and geophysical data or regional structure numerical simulations are required to determine how the subduction of the early Pacific Plate or the Australian Plate affects the deep mantle materials of the Bismarck Basin, as well as to investigate the relationship between the dynamic factors driving the expansion of the Manus Basin.

Acknowledgements

We thank Wessel & Smith (https://www.soest.hawaii.edu/ gmt/) for the free use of GMT software to produce most figures in this work. We are grateful for the codes of double-difference tomography shared by Prof. Haijiang Zhangat the University of Science and Technology of China. This work was supported by the National Natural Science Foun- dation of China (Nos. 41906048 and 91858215).

Abers, G. A., and Roecker, S. W., 1991. Deep structure of an arc-continent collision: Earthquake relocation and inversion for up- per mantle P and S wave velocities beneath Papua New Guinea., 96 (B4): 6379-6401.

Binns, R. A., and Scott, S. D., 1993. Activity forming polymetallic sulfide deposits associated with felsic volcanic rocks in the eastern Manus Back-arc basin, Papua New Guinea., 88 (8): 2226-2236.

Both, R., Crook, K., Taylor, B., Brogan, S., Chappell, B., Frankel, E.,., 1986. Hydrothermal chimneys and associated fauna in the Manus Back-arc basin, Papua New Guinea., 67 (21): 489-490.

DeMets, C., Gordon, R. G., and Argus, D. F., 2010. Geologically current plate motions., 181 (1): 1-80.

Feng, B., Zhang, H., Gong, W., Xing, J., Li, D., Xu, C.,., 2022. Teleseismic P-wave tomography of the New Guinea- So- lomon arc system., 21 (3): 694-706.

Finlayson, D. M., Cull, J. P., Wiebenga, W. A., Furumoto, A. S., and Webb, J. P., 1972. New Britain-New Ireland crustal seismic refraction investigations 1967 and 1969., 29 (3): 245-253.

Franz, L., and Romer, R. L., 2010. Different styles of metasomatic veining in ultramafic xenoliths from the Tubaf seamount (Bis- marck microplate, Papua New Guinea)., 114: 30-53.

Gong, W., Jiang, X., Xing, J., Li, D., and Xu, C., 2019. Subduction dynamics of the New-Guinea-Solomon arc system: Constraints from the subduction initiation of the plate., 39 (9): 115-129 (in Chinese with English abstract).

Heidbach, O., Tingay, M., Barth, A., Reinecker, J., Kurfe, D., and Müller, B., 2010. Global crustal stress pattern based on the World Stress Map database release 2008., 482 (1): 3-15.

Holm, R. J., and Richards, S. W., 2013. A re-evaluation of arc- continent collision and along-arc variation in the Bismarck Sea region, Papua New Guinea., 60 (5): 605-619.

Holm, R. J., Rosenbaum, G., and Richards, S. W., 2016. Post 8 Ma reconstruction of Papua New Guinea and Solomon Islands: Microplate tectonics in a convergent plate boundary setting., 156: 66-81.

Honza, E., Miyazaki, T., and Lock, J., 1989. Subduction erosion and accretion in the Solomon Sea region., 160 (1-4): 49-62.

Hyndman, R. D., and Peacock, S. M., 2003. Serpentinization of the forearc mantle., 212 (3-4): 417-432.

Johnson, T., and Molnar, P., 1972. Focal mechanisms and plate tectonics of the Southwest Pacific., 77 (826): 5000-5032.

Lee, S. M., and Ruellan, E., 2006. Tectonic and magmatic evolution of the Bismarck Sea, Papua New Guinea: Review and new synthesis. In:. Chirstie, D. M.,., eds., American Geophysical Union, Washington, D. C., 263-286.

Ma, Y., Yin, X., Wang, X., Chen, S., and Zeng, Z., 2018. Influence of tectonic-magmatism on hydrothermal activity: A case study of the Manus Basin., 42 (5): 163-171 (in Chinese with English abstract).

Ma, Y., Zeng, Z. G., Chen, S., Yin, X., and Wang, X., 2017. Origin of the volcanic rocks erupted in the Eastern Manus Basin: Basaltic andesite-andesite-dacite associations., 16 (3): 389-402.

McGuire, J. D., and Wiens, D. A., 2013. A double seismic zone in New Britain and the morphology of the Solomon Plate at intermediate depths., 22 (15): 1965-1968.

Mori, J., 1989. The New Ireland earthquake of July 3, 1985 and associated seismicity near the Pacific Solomon Sea-Bismarck Sea triple junction., 55 (1-2): 144-153.

Mu, D., Li, S., Suo, Y., Zhu, J., and Yu, S., 2019. Tectonic and geodynamic mechanism of back-arc-rifting derived micro- blocks: Insights from back-arc spreading in the West Pacific., 43 (4): 665-677 (in Chinese with English abstract).

Ortega-Osorio, A., and Scott, S. D., 2001. Morphological and che-mical characterization of neutrally buoyant plume-derived particles at the eastern Manus Basin hydrothermal field, Papua New Guinea., 39: 17-31.

Paige, C. C., and Saunders, M. A., 1982. LSQR: An algorithm for sparse linear equations and sparse least squares., 8: 43-71.

Pegler, G., Das, S., and Woodhouse, J. H., 1995. A seismological study of the eastern New Guinea and the western Solomon Searegions and its tectonic implications., 122 (3): 961-981.

Sinton, J. M., Ford, L. L., Bruce, C., and Mcculloch, M. T., 2003. Magma genesis and mantle heterogeneity in the manus back- arc basin, Papua New Guinea., 44 (1): 159- 195.

Spakman, W., Lee, S., and Hilst, R., 1993. Travel-time tomography of the European-Mediterranean mantle down to 1400km., 79 (1-2): 3-74.

Stern, R. J., 2002. Subduction zones., 40 (4): 3-38.

Tatsumi, Y., and Eggins, S., 1995.. Blackwell Science, Boston, 211pp.

Tatsumi, Y., Hamilton, D. L., and Nesbitt, R. W., 1986. Chemical characteristics of fluid phase released from a subducted li- thosphere and origin of arc magmas: Evidence from high-pres- sure experiments and natural rocks., 29: 293-309.

Tatsumi, Y., Sakuyama, M., Fukuyama, H., and Kushiro, I., 1983. Generation of arc basalt magmas and thermal structure of the mantle wedge in subduction zones (Japan arc)., 88: 5815-5825.

Taylor, B., 1979. Bismarck Sea–Evolution of a back-arc basin., 7 (4): 171-174.

Thal, J., Tivey, M., Yoerger, D., Jöns, N., and Bach, W., 2014. Geo- logic setting of PACManus hydrothermal area–High resolu- tion mapping andobservations., 355: 98-114.

Tregoning, P., Jackson, R. J., McQueen, H., Lambeck, K., Stevens, C., Little, R. P.,., 1999. Motion of the South Bismarck Plate, Papua New Guinea., 26 (23): 3517-3520.

Tregoning, P., Lambeck, K., Stolz, A., Morgan, P., McClusky, S. C., van ber Beek, P.,., 1998. Estimation of current plate motions in Papua New Guinea from global positioning system observation., 103 (B6): 12181-12203.

Waldhauser, F., and Ellsworth, W., 2000. A double-difference earthquake location algorithm: Method and application to the northern hayward fault California., 90: 1353-1368.

Wallace, L. M., Stevens, C., Silver, E., McCaffrey, R., Loratung, W., Hasiata, S.,, 2004. GPS and seismological constraints on active tectonics and arc-continent collision in Papua New Guinea: Implications for mechanics of microplate rotations in a plate boundary zone., 109 (B5): B05404.

Weiler, P. D., and Coe, R. S., 2000. Rotations in the actively colliding Finisterre arc terrane: Paleomagnetic constraints on plio-Pleistocene evolution of the South Bismarck microplate, northeastern Papua New Guinea., 316: 297-325.

Whitmore, G. P., Johnson, D. P., Crook, K. A. W., Galewsky, J., and Silver, E. A., 1997. Convergent margin extension associatedwith arc-continent collision: The Finsch Deep, Papua New Gui- nea., 16 (1): 77-87.

Wyllie, P. J., 1988. Magma genesis, plate tectonics, and chemical differentiation of the earth., 26 (3): 370-404.

Wyss, M., Hasegawa, A., and Nakajima, J., 2001. Source and path of magma for volcanoes in the subduction zone of northeastern Japan., 28 (9): 1819-1822.

Yamasaki, T., and Seno, T., 2003. Double seismic zone and dehydration embrittlement of the subducting slab., 108 (B4): 2212, DOI: 10.10 29/2002JB001918.

Zhang, H. J., and Thurber, C., 2003. Double-difference tomography: The method and its application to the Hayward fault, California., 93 (5): 1875-1889.

Zhang, H. J., and Thurber, C., 2006. Development and applications of double-difference seismic tomography., 163 (2-3): 373-403.

Zhang, K. L., and Wei, D. P., 2008. Progresses of the researches and the causing mechanisms on the double seismic zones with- in the subduction zones around the Pacific Ocean., 23 (1): 31-39 (in Chinese with English abstract).

Zhao, H., Zeng, Z., Yin, X., and Chen, S., 2014. Silicon and oxygen isotopic composition of igneous rocks from the eastern Ma- nus Basin., 13 (3): 421-427.

Zoback, M. L., 1992. First- and second-order patterns of stress in the lithosphere: The world stress map project., 97 (B8): 11703-11728.

(December 3, 2021; revised May 20, 2022; accepted July 22, 2022)

© Ocean University of China, Science Press and Springer-Verlag GmbH Germany 2023

Corresponding author. E-mail: lichaoyang@sdust.edu.cn

(Edited by Chen Wenwen)