William J Shinevar

Hannah F Mark

Fiona Clerc

Emmanuel A Codillo

Jianhua Gong

Jean-Arthur Olive

Stephanie M Brown

Paris T Smalls

Yang Liao

Véronique Le Roux

Mark D Behn

Causes of oceanic crustal thickness oscillations along a 74-Myr Mid-Atlantic Ridge flow line

We conducted one of the longest continuous geophysical surveys along a 74-Myr spreading-parallel flow line across the Mid-Atlantic Ridge.

2. Spectral densities of bathymetry and gravity data show concurrent peaks at 390-, 550-, and 950-kyr periods, and diffuse power at >1 Myr.

3.

A negative correlation between fault spacing and gravity-derived crustal thickness suggests a link between magma input and fault style.

Introduction

The oceanic lithosphere covers more than two-thirds of the Earth's surface and plays an essential role in plate tectonics. Oceanic lithosphere forms at mid-ocean ridges (MORs) through a combination of magmatic and tectonic activity [START_REF] Macdonald | Volcanic growth faults and the origin of Pacific abyssal hills[END_REF][START_REF] Searle | Mid-ocean ridges[END_REF]. Variability in the magma supply associated with differences in spreading rate, axial thermal structure, and/or mantle melting processes has been linked to first-order differences in ridge morphology (e.g., Chen & Morgan, 1990a;[START_REF] Chen | A nonlinear rheology model for mid-ocean ridge axis topography[END_REF]Ito & Behn, 2008;[START_REF] Roth | Discrete characterization of abyssal hills: Unraveling temporal variations in crustal accretion processes at the 10° 30′ N segment, East Pacific Rise[END_REF], fault behavior (e.g., [START_REF] Buck | Modes of faulting at mid-ocean ridges[END_REF]Behn & Ito, 2008;[START_REF] Tucholke | Role of melt supply in oceanic detachment faulting and formation of megamullions[END_REF], and the chemistry of mid-ocean ridge basalts (e.g., [START_REF] Bonatti | Mantle thermal pulses below the Mid-Atlantic Ridge and temporal variations in the formation of oceanic lithosphere[END_REF][START_REF] Behn | Melting systematics in mid-ocean ridge basalts: Application of a plagioclase-spinel melting model to global variations in major element chemistry and crustal thickness[END_REF].

Temporal variability in magmatic and tectonic forcing at mid-ocean ridges has been hypothesized and/or observed on timescales from millions to thousands of years.

For example, [START_REF] Bonatti | Mantle thermal pulses below the Mid-Atlantic Ridge and temporal variations in the formation of oceanic lithosphere[END_REF] found a long-term increase in gravity-derived crustal thickness along the Vema Transform Fault (Guinea Fracture Zone) in the Atlantic, which they hypothesized to result from a long-term (>20 Myr) increase in mantle potential temperature. Previous studies in the Mid-Atlantic also found 2-4 Myr period variations in gravity-derived crustal thickness, which were interpreted as representing changes in mantle upwelling and flow [START_REF] Bonatti | Mantle thermal pulses below the Mid-Atlantic Ridge and temporal variations in the formation of oceanic lithosphere[END_REF][START_REF] Pariso | Temporal and spatial variations in crustal accretion along the Mid-Atlantic Ridge (29°-31° 30′ N) over the last 10 my: Implications from a three-dimensional gravity study[END_REF] or as reflecting alternating phases of dominantly magmatic and amagmatic extension [START_REF] Tucholke | Segmentation and crustal structure of the western Mid-Atlantic Ridge flank, 25° 25′-27° 10′ N and 0-29 my[END_REF]. Further, [START_REF] Minshull | Crustal structure of the Southwest Indian Ridge at 66 E: Seismic constraints[END_REF] attributed seismically derived crustal thickness variations with periods of ~3 Myr near the ultraslow-spreading Southwest Indian Ridge to be the result of episodic melt flow and tectonic extension during amagmatic periods.

On slightly shorter time-scales, Canales et al. (2000a) reported seismicallyderived crustal thickness variations with periods of ~0.5-1 Myr near the axis of the Mid-Atlantic Ridge close to the Kane fracture zone, which they attributed to temporal variability in magma supply. Measured cross sections of the Oman ophiolite also depict thickness variations over 5-10 km scale attributed to accretion above mantle diapirs [START_REF] Nicolas | Variable crustal thickness in the Oman ophiolite: implication for oceanic crust[END_REF]. Assuming a half spreading rate between 25 and 50 mm yr -1 for the crust formed in the Oman ophiolite, the 5-10 km scale thickness variations correspond to variations with a periodicity of 100-400 kyr. High-resolution sampling across the intermediate-spreading Juan de Fuca Ridge [START_REF] Ferguson | A 65 ky time series from sediment-hosted glasses reveals rapid transitions in ocean ridge magmas[END_REF] revealed a rapid change (<10 kyr) in basalt geochemistry coincident with 1-km increase in crustal thickness [START_REF] Carbotte | Variable crustal structure along the Juan de Fuca Ridge: Influence of on-axis hot spots and absolute plate motions[END_REF], demonstrating that changes in magma supply can also occur on very short time scales. Some short period (<100 kyr) variations in crustal production have been proposed

to result from external forcing on the MOR magmatic system. [START_REF] Crowley | Glacial cycles drive variations in the production of oceanic crust[END_REF] hypothesized that fluctuations in sea level caused by Milanković cycles (23, 41, and 100 kyr) could modulate the amount of magma produced by adiabatic decompression melting below MORs by modulating mantle pressure. [START_REF] Tolstoy | Mid-ocean ridge eruptions as a climate valve[END_REF] also suggested that orbital forcing affects crustal stresses and thus the dynamics of melt extraction on Milanković time scales, particularly at fast-spreading ridges.

One challenge in deciphering temporal changes in the structure and evolution of oceanic lithosphere is that most previous studies have focused on the formation of young lithosphere at or near the ridge axis. For example, to understand the magmatic input to mid-ocean ridges, studies have focused on variations in the chemistry of mid-ocean ridge basalts (e.g., [START_REF] Langmuir | Petrological systematics of midocean ridge basalts: Constraints on melt generation beneath ocean ridges. Mantle flow and melt generation at mid-ocean ridges[END_REF][START_REF] Gale | The mean composition of ocean ridge basalts[END_REF], near-axis seismic and gravity measurements of crustal thickness (e.g., [START_REF] Kuo | Gravity anomalies of the ridge-transform system in the South Atlantic between 31 and 34.5 S: Upwelling centers and variations in crustal thickness[END_REF][START_REF] Detrick | Gravity anomalies and crustal thickness variations along the Mid-Atlantic Ridge between 33° N and 40° N[END_REF][START_REF] Hooft | Crustal thickness and structure along three contrasting spreading segments of the Mid-Atlantic Ridge, 33.5°-35[END_REF], and the physical dynamics of magma chambers and their eruptions (e.g., [START_REF] Sinton | Mid-ocean ridge magma chambers[END_REF][START_REF] Rubin | Dike-induced faulting in rift zones of Iceland and Afar[END_REF][START_REF] Tan | Dynamics of a seafloor-spreading episode at the East Pacific Rise[END_REF]. Near-axis geophysical observables, such as gravity, bathymetry, and fault style, have also been analyzed to understand the way magmatism and tectonics interact to shape the ridge axis (e.g., [START_REF] Kuo | Gravity anomalies of the ridge-transform system in the South Atlantic between 31 and 34.5 S: Upwelling centers and variations in crustal thickness[END_REF][START_REF] Detrick | Gravity anomalies and crustal thickness variations along the Mid-Atlantic Ridge between 33° N and 40° N[END_REF][START_REF] Hooft | Crustal thickness and structure along three contrasting spreading segments of the Mid-Atlantic Ridge, 33.5°-35[END_REF][START_REF] Buck | Modes of faulting at mid-ocean ridges[END_REF][START_REF] Howell | Magmatic and tectonic extension at the Chile Ridge: Evidence for mantle controls on ridge segmentation[END_REF]2019). By contrast, relatively few studies have investigated off-axis oceanic lithosphere (e.g., [START_REF] Morris | Seismic structure of oceanic crust in the western North Atlantic[END_REF][START_REF] Tucholke | Segmentation and crustal structure of the western Mid-Atlantic Ridge flank, 25° 25′-27° 10′ N and 0-29 my[END_REF]Lizzaralde et al., 2004) where mass wasting and sedimentation restrict sampling and complicate the bathymetric record of faulting. Thus, temporal variability has typically been inferred from very short acrossaxis studies (e.g., Canales et al., 2000a;[START_REF] Crowley | Glacial cycles drive variations in the production of oceanic crust[END_REF], as well as comparison of spatial variations between different ridge segments of similar spreading rates (e.g. [START_REF] Purdy | Relationship between spreading rate and the seismic structure of mid-ocean ridges[END_REF].

To better quantify the temporal variability in magma supply and fault style in oceanic lithosphere over timescales of 10 4 -10 7 years, we acquired a 1400-km-long geophysical transect collecting bathymetry, gravity, and magnetic data during the SCARF cruise (Student-led Cruise Along a Ridge Flow Line) (#AR023-02) aboard the R/V Neil Armstrong in October 2017. The cruise track followed a single, spreading-parallel flow line between 30˚ and 45˚W crossing the Mid-Atlantic Ridge (MAR) at 35.8˚N, ~100 km north of the Oceanographer Fracture Zone (Figure 1). The cruise yielded one of the longest geophysical surveys along a single ridge flow line and covered ~28 Ma lithosphere on the African plate to 74 Ma lithosphere on the North American plate. Our analysis of this flow line allows us to estimate key properties of the oceanic lithosphere, including crustal thickness, effective elastic plate thickness, and fault characteristics to provide constraints on the tectonic and magmatic processes that control temporal variation of crustal formation at the same mid-ocean ridge segment.

Geologic Setting

The Central Atlantic Ocean, the portion of the Atlantic Ocean south of the Azores plateau and north of Guinean Fracture Zone, began spreading at ~190 Ma [START_REF] Labails | An alternative early opening scenario for the Central Atlantic Ocean[END_REF], likely triggered by the 203 and 193 Ma surges of Central Atlantic magmatic province volcanism [START_REF] Jourdan | 40Ar/39Ar ages of CAMP in North America: implications for the Triassic-Jurassic boundary and the 40K decay constant bias[END_REF]. Early accretion in the Central Atlantic Ocean was highly asymmetric, adding significantly more material to the North American plate.

The Central Atlantic Ocean originally accreted slowly (half spreading rate of 4 mm yr -1 ) until 170 Ma, when the spreading direction shifted from northwest-southeast to east-west and spreading rate doubled. Afterwards, the half spreading rate remained around 10 mm yr -1 with periods of faster rates up to 30 mm yr -1 [START_REF] Labails | An alternative early opening scenario for the Central Atlantic Ocean[END_REF][START_REF] Seton | Global continental and ocean basin reconstructions since 200 Ma[END_REF].

The Mid-Atlantic Ridge north of the Oceanographer Fracture Zone and south of the Tydeman Fracture Zone consists of 20-45 km long segments separated by nontransform offsets. The SCARF cruise track followed a flow line that crossed the ridge axis near the center of the third segment north of the Oceanographer Fracture Zone at 35.8˚N, 34.2˚W. The crossing point lies approximately 90 km north of the Oceanographer Fracture Zone and 260 km south of Tydeman Fracture Zone. This segment has been named previously as the South Alvin Mid-Atlantic Ridge (S. AMAR) (e.g., [START_REF] Detrick | Gravity anomalies and crustal thickness variations along the Mid-Atlantic Ridge between 33° N and 40° N[END_REF] and AMARR91 by [START_REF] Gale | The mean composition of ocean ridge basalts[END_REF]. The S. AMAR segment is 40-km long and has a continuous, hourglass-shaped axial rift valley [START_REF] Paulatto | Heterogeneous and asymmetric crustal accretion: New constraints from multibeam bathymetry and potential field data from the Rainbow area of the Mid-Atlantic Ridge (36˚ 15'N)[END_REF][START_REF] Eason | Segment-scale variations in seafloor volcanic and tectonic processes from multibeam sonar imaging, Mid-Atlantic Ridge Rainbow region (35˚ 45'-36˚ 35'N)[END_REF].

The surrounding near-axis region of the Mid-Atlantic Ridge has been well studied. The Rainbow Hydrothermal Vent Field [START_REF] German | The Rainbow hydrothermal plume, 36˚15[END_REF] is situated approximately 50 km northeast of the S. AMAR segment, and the French-American Mid-Ocean Undersea Study (Project FAMOUS) [START_REF] Ballard | Manned submersible observations in the FAMOUS area: Mid-Atlantic Ridge[END_REF] investigation area was located two segments further to the north. The axial rift valley is shallow where it intersects the SCARF cruise track, approximately 2580 m deep and 8 km wide. The segment currently has a slightly asymmetric ridge-perpendicular half spreading rate of 14 mm yr -1 eastward and 13 mm yr -1 westward as calculated from magnetic isochrones in this study (see Section 3.1). Rare-earth element and isotopic chemistry of dredged MORB samples for this segment, as well as other nearby segments show heterogeneity, suggesting at least three mantle components, one of which comes from the Azores plume [START_REF] White | The nature and origin of geochemical variation in Mid-Atlantic Ridge basalts from the central North Atlantic[END_REF][START_REF] Bougault | Mantle heterogeneity from trace elements: MAR triple junction near 14 N[END_REF][START_REF] Shirey | Three-component isotopic heterogeneity near the Oceanographer transform, Mid-Atlantic Ridge[END_REF].

Observations: Underway Geophysical Data Collection

Magnetics

We used sea surface total magnetic field data to constrain seafloor age along the SCARF cruise track. The data were acquired using a SeaSPY Overhauser magnetometer towed at an offset of 300 meters behind the ship. The raw sea surface total magnetic field data was corrected for (1) navigation offset from ship-to-magnetometer layback, (2) outlying data points, (3) the international geomagnetic reference field (IGRF-11 model, [START_REF] Finlay | International geomagnetic reference field: the eleventh generation[END_REF], and (4) diurnal field variations (relative to Fredericksburg, Virginia, USA, (FRD)). Sea surface crustal magnetic anomalies were obtained by correcting the total field data for the paleo-inclinations (30.5°-35.2°) and -declinations (49.7°-54.63°) estimated for the Atlantic using the method of [START_REF] Schettino | Quantitative plate tectonics[END_REF] and the paleopoles of [START_REF] Schettino | Apparent polar wander paths for the major continents (200 Ma to the present day): a palaeomagnetic reference frame for global plate tectonic reconstructions[END_REF]. The resulting magnetic anomalies along the SCARF transect are in good agreement with the EMAG2v3 (version 3) global magnetic anomaly grid [START_REF] Meyer | EMAG2: Earth magnetic anomaly grid (2arc-minute resolution) version 3[END_REF] when smoothed (Figure A1) and with magnetic anomaly measurements from nearby or crossing ship paths from NCEI's Marine Trackline Geophysical database (https://www.ngdc.noaa.gov/mgg/geodas/trackline.html) (Figure 1b).

We estimated the ages represented by the magnetic anomalies along the SCARF transect by matching the observed anomalies to a calculated synthetic magnetic anomaly profile. The synthetic model assumes a magnetization square wave based on a magnetized layer at a constant depth of 3.5 km. The layer is 1 km thick and has a magnetization of ±5 A/m with polarities based on the [START_REF] Müller | Age, spreading rates, and spreading asymmetry of the world's ocean crust[END_REF] ocean age database and the [START_REF] Wei | Revised age calibration points for the geomagnetic polarity time scale[END_REF] Geomagnetic Polarity Time Scale. We filtered this synthetic magnetization profile with a 2-km Gaussian window to account for emplacement effects, following the approach of [START_REF] Tivey | Magnetization of 0-29 Ma ocean crust on the Mid-Atlantic Ridge, 25° 30′ to 27° 10′[END_REF]. The synthetic magnetization was converted to a synthetic magnetic anomaly profile using the Parker (1972) Fourier summation approach.

To assign isochron numbers and ages to magnetic anomalies along the SCARF cruise track, we visually correlated the synthetic magnetic anomaly profile with the measured magnetic anomalies (Figure 2). We defined the ridge axis as the midpoint of the axial valley and set its age to 0 Ma. Spreading rates were predicted on the basis of the correlated magnetic anomalies and agree well with those calculated along the SCARF transect using the global seafloor age model of [START_REF] Müller | Age, spreading rates, and spreading asymmetry of the world's ocean crust[END_REF] (Figure 3d).

Bathymetry and Fault Identification

Bathymetry data (Figure 3a) were collected using a Kongsberg EM122 12-kHz multibeam echosounder. This multibeam echosounder has 0.2% error directly below the ship with a maximal 0.6% error at 60-70˚. Center beam bathymetry was collected at 1 Hz and multi-beam bathymetry was collected at 0.33 Hz. The data were cleaned for outlying data points and gridded with the MB-System software [START_REF] Caress | MB-System: Mapping the Seafloor[END_REF]. Figure A2.1-18 shows images of all the trackline multi-beam bathymetry (bathymetry data available through Marine Geoscience Data System). The swath width is ~6 times the water depth. We identified individual faults manually using a combination of the center beam and multibeam bathymetry. We identified axis-facing fault scarps as changes in topographic gradient observed in the center beam bathymetry data that correspond to linear features in the multibeam bathymetry. Abyssal hills may not be entirely faultcontrolled, and some may be relict volcanic ridges. Volcanic ridges should manifest as symmetric ridges, displaying both inward-and outward-facing slopes within a few km of one another. Although common at fast spreading rates, volcanic ridges are typically not preserved at slower spreading rates [START_REF] Carbotte | Comparison of seafloor tectonic fabric at intermediate, fast, and super fast spreading ridges: Influence of spreading rate, plate motions, and ridge segmentation on fault patterns[END_REF], and such morphological features were not seen along the transect. Here we interpret all fault scarp picks as abyssal hill-bounding faults. We estimate the throw (vertical offset) and heave (horizontal offset) of each fault by picking the shallowest depth of the fault in the center beam bathymetry as the scarp top and the deepest depth of the fault as the scarp bottom.

We only picked faults that have an identifiable throw of greater than 10 m. Using this approach, we identified 415 fault scarps along the SCARF transect (red dots in Figure 3a). Individual fault statistics are available in Supplementary Table 1. We find most fault slopes have less than 30˚ angles, likely due to gravitational mass wasting of the fault scarps (Cannat et al., 2013) (Figure A3). We calculated fault spacing (bin length divided by the number of picked faults in a bin, Figure 3b) and total fault heave in overlapping 100 km bins every 25 km along the flow line. There are 12 or more faults in each bin.

The resulting bin spacing precludes analysis of short-period (< 2 Myr) variations in fault spacing.

We validated our fault picks by examining the frequency distribution of fault throws (Figure 4). Natural fault populations are predicted to display an exponential frequency distribution for faults having lengths comparable to or greater than the brittle layer thickness [START_REF] Cowie | Fault strain and seismic coupling on mid-ocean ridges[END_REF][START_REF] Carbotte | Comparison of seafloor tectonic fabric at intermediate, fast, and super fast spreading ridges: Influence of spreading rate, plate motions, and ridge segmentation on fault patterns[END_REF] and exponential distributions in scarp heights and fault spacing have been observed at the East Pacific Rise and Chile Ridge [START_REF] Bohnenstiehl | Faulting and fault scaling on the median valley floor of the trans-Atlantic geotraverse (TAG) segment, 268N on the Mid-Atlantic Ridge[END_REF]Bohnsenstiehl & Carbotte, 2001;[START_REF] Howell | Magmatic and tectonic extension at the Chile Ridge: Evidence for mantle controls on ridge segmentation[END_REF]. Consistent with these results, we observe that scarp frequency decays exponentially with scarp throw for our picks, suggesting that our scarp database is sampling the true fault population. The exponential decay of fault population fades at large throw (>0.7 km); we only find three faults with throw greater than 0.7 km.

A major difference between previous fault population studies and ours is that most previous studies analyzed only near axis faults along multiple ridge segments, while most of our picked faults are located far from the ridge axis in >10 Ma lithosphere. For instance, the furthest off axis fault analyzed by [START_REF] Howell | Magmatic and tectonic extension at the Chile Ridge: Evidence for mantle controls on ridge segmentation[END_REF] is 100 km away from the ridge axis, corresponding to ~3 Ma lithosphere. To check the validity of our off-axis fault picks, we compared the distributions of fault throw for (1) all faults within 200 km of the ridge, (2) all faults within 400 km of the ridge, and (3) the entire transect. We find no change of slope in the cumulative throw distributions between those calculated for 200 and 400 km cut-offs below a throw of ~0.6 km. The change in slope for throw > 0.6 km is likely due to the low number of high throw faults found within 400 km of the ridge axis along our transect, making it difficult for a meaningful comparison. Beyond 400 km offaxis, we identify no faults with throw > 0.5 km. In terms of total fault occurrence, we find 138 faults (33% of the total faults found along the transect) at distances >400 km offaxis (45% of our track line). The lack of identified faults off-axis is likely due to sedimentation and mass wasting reworking seafloor topography. This interpretation is supported by the observation that at higher sediment thicknesses (>400 m) there is an apparent increase in fault spacing, while for lower sediment thicknesses (< ~300 m) nearer to the ridge axis (<400 km away), fault spacing does not correlate with sediment thickness (Figure 3b). Thus, we limit our analyses of fault statistics on this transect to within 400 km of the ridge and interpret fault spacing variations within this distance to reflect changes in tectonic behavior.

We calculated the fraction of spreading accommodated by magmatic accretion, M (0 for purely tectonic spreading and 1 for purely magmatic spreading; [START_REF] Buck | Modes of faulting at mid-ocean ridges[END_REF] in each 100-km bin by subtracting the cumulative heave from the bin length and dividing by the bin length. M values calculated in this manner are maxima, as apparent fault heave recorded by scarp morphology is reduced through sedimentation and mass wasting.

In general, these effects should become more significant in older lithosphere farther offaxis. Focusing on bins within 400 km of the ridge axis, where we have the best estimates of fault statistics, we find that M varies from 0.6-0.85 and is smallest near the ridge (Figure 3c). Farther off-axis, variations in fault spacing are likely to be better proxies for variations in tectonic behavior, as spacing will only be affected by mass wasting and sedimentation if an existing fault is degraded or buried such that it can no longer be identified in the center beam data.

Gravity

Shipboard gravity data was collected at a 1 Hz sampling rate along the SCARF transect using a Bell BGM-3 gravimeter. Corrections were applied for (1) instrument calibration, (2) instrument drift, (3) latitude variations, and (4) the Eötvös effect. We filtered the corrected gravity data with a 2-pole low-pass Butterworth filter with corner frequency 0.004 per sample and removed outlying data points. The free-air anomaly (FAA) was calculated by subtracting the reference gravity along the GRS80 ellipsoid [START_REF] Moritz | Geodetic reference system 1980[END_REF] from the observed gravity. Our data is in excellent agreement with the 1arcminute satellite-derived marine FAA version 24 [START_REF] Sandwell | New global marine gravity model from CryoSat-2 and Jason-1 reveals buried tectonic structure[END_REF] (Figure 1c).

The mantle Bouguer anomaly (MBA) was calculated by subtracting the water-sediment, sediment-crust, and crust-mantle density interfaces from the FAA following the spectral method of [START_REF] Parker | The rapid calculation of potential anomalies[END_REF]. We assume a 6-km thick crustal layer. Sediment thickness (Figure 3b) was taken from [START_REF] Divins | Total sediment thickness of the world's oceans & marginal seas[END_REF], which interpolated sediment thickness from compiled isopach maps, ocean drilling records, and seismic reflection data. To reduce edge effects in the calculation of the MBA, we added 500 km of satellite-derived data [START_REF] Sandwell | New global marine gravity model from CryoSat-2 and Jason-1 reveals buried tectonic structure[END_REF] to each end of the SCARF transect. Densities used to calculate density contrasts for water, sediment, crust, and mantle were assumed to be 1000, 2100, 2800, and 3300 kg m -3 , respectively, with the sediment density taken from [START_REF] Hamilton | Sound velocity-density relations in sea-floor sediments and rocks[END_REF]. To obtain the residual mantle Bouguer anomaly (RMBA) (Figure 3e), we further corrected the MBA for thermal contraction of the oceanic lithosphere assuming a halfspace cooling model [START_REF] Turcotte | Geodynamics[END_REF] with the magnetically derived ages (Figure 2), 𝑇 𝑚𝑎𝑛𝑡𝑙𝑒 = 1350 ˚C, 𝑇 𝑠𝑢𝑟𝑓𝑎𝑐𝑒 = 0 ˚C, thermal expansivity 𝛼 = 5x10 -5 ˚C-1 , and thermal diffusivity 𝜅 = 10 -6 m 2 s -1 .

Gravity-Derived Crustal Thickness

The RMBA can be used to estimate relative variations in crustal thickness assuming a constant crustal density. Relative variations in crustal thickness derived from this approach have been shown to compare well to constraints from co-located seismic studies (e.g., [START_REF] Tolstoy | Crustal thickness on the Mid-Atlantic Ridge: Bull's-eye gravity anomalies and focused accretion[END_REF][START_REF] Hooft | Crustal thickness and structure along three contrasting spreading segments of the Mid-Atlantic Ridge, 33.5°-35[END_REF][START_REF] Wang | Crustal thickness anomalies in the North Atlantic Ocean basin from gravity analysis[END_REF]. Here we calculated relative crustal thickness variations by downward continuing the de-meaned RMBA to a constant depth of 10 km below the sea surface [START_REF] Kuo | Gravity anomalies of the ridge-transform system in the South Atlantic between 31 and 34.5 S: Upwelling centers and variations in crustal thickness[END_REF][START_REF] Wang | Crustal thickness anomalies in the North Atlantic Ocean basin from gravity analysis[END_REF]. To remove short-wavelength noise, we filtered the RMBA before downward continuation, using a cosine taper designed such that signals with wavelengths less than 25 km were filtered out and signals with wavelengths greater than 45 km were fully retained. Because the average half spreading rate along the SCARF transect is ~14 km Myr -1 , this removes crustal thickness variations occurring at wavelengths less than ~2 Myr. We then calculated crustal thickness variations with respect to a 6 km crust (Figure 3f), corresponding to the average crustal thickness obtained at a nearby ridge segment using seismic refraction data [START_REF] Hooft | Crustal thickness and structure along three contrasting spreading segments of the Mid-Atlantic Ridge, 33.5°-35[END_REF]. Our results agree well with the regional estimates of gravity-derived crustal thickness calculated by [START_REF] Wang | Crustal thickness anomalies in the North Atlantic Ocean basin from gravity analysis[END_REF].

Kafka Dome

Along our transect, we observed a dome-shaped bathymetric high at 35˚ 34' N, 32˚ 54' W (Figure 5), which we name Kafka Dome. Kafka Dome was identified as an oceanic core complex based on its dome morphology with a topographic prominence of ~1 km and a ~7 km length widening away from the axis into an elongated ridge. This morphology is consistent with the breakaway ridge geometry observed at well-developed oceanic core complexes at 13˚20'N and 13˚30'N on the MAR (Figure 5c, [START_REF] Escartín | Tectonic structure, evolution, and the nature of oceanic core complexes and their detachment fault zones (13 20' N and 13 30' N, Mid Atlantic Ridge)[END_REF] and Mt. Dent on the Mid-Cayman spreading center [START_REF] Hayman | Oceanic core complex development at the ultraslow spreading Mid-Cayman Spreading Center[END_REF]. The measured size and topographic height of the dome are likely minimal estimates of the original dome structure due to off-axis sedimentation and mass wasting. We also observe a local ~8 mGal increase in the RMBA skewed towards the ridge axis (Figure 5b).

Similar RMBA increases at other oceanic core complexes are commonly interpreted as positive mass anomalies caused by thinner crust in the fault footwall [START_REF] Tucholke | Megamullions and mullion structure defining oceanic metamorphic core complexes on the Mid-Atlantic Ridge[END_REF]2008). We do not observe corrugations as seen at many near ridge core complexes (e.g., [START_REF] Cann | Corrugated slip surfaces formed at ridge-transform intersections on the Mid-Atlantic Ridge[END_REF], likely due to off-axis sedimentation. Taken together, these observations point toward Kafka Dome as a rare example of an inactive off-axis oceanic core complex.

The discovery of a fossil core complex indicates that within the last 10 Ma, faulting style at the S. AMAR segment has changed from an oceanic core complex regime to a regularly spaced abyssal hill regime characteristic of presently accreting seafloor in the vicinity of the Oceanographer Fracture Zone [START_REF] Escartín | Central role of detachment faults in accretion of slow-spreading oceanic lithosphere[END_REF]. The Rainbow Massif and Pot of Gold Massif are the closest known large-offset detachment faults [START_REF] Paulatto | Heterogeneous and asymmetric crustal accretion: New constraints from multibeam bathymetry and potential field data from the Rainbow area of the Mid-Atlantic Ridge (36˚ 15'N)[END_REF][START_REF] Eason | Segment-scale variations in seafloor volcanic and tectonic processes from multibeam sonar imaging, Mid-Atlantic Ridge Rainbow region (35˚ 45'-36˚ 35'N)[END_REF]. The coincidence of Kafka Dome with the region of larger fault spacing located at 8-22 Myr on both sides of the ridge axis (shaded bars in Figure 3) suggests that this interval may be characteristic of a long-lived oceanic core complex regime.

Time-Series Analysis

Using the magnetically derived ages and our shipboard geophysical datasets, we apply spectral analysis to decompose the bathymetry and gravity data in space and time and determine what frequency range(s) hold the most power in our data. Because various tectonic and magmatic processes act on different timescales, the spectral power density can distinguish the relative importance of these processes. We also calculate the coherence and admittance between the gravity and simultaneously-observed center beam bathymetry data, which we use to estimate the elastic plate thickness, 𝑇 𝑒 , along the flow line.

Spectral Analysis

Time-domain spectral analysis on the bathymetry and RMBA (Figure 6) was performed using the multitaper method [START_REF] Percival | Spectral Analysis for Physical Applications: Multitaper and Conventional Univariate Techniques[END_REF] with seven tapers.

We used ages defined by linear interpolation of our magnetic isochron picks. We investigate spectral power in the time domain such that spreading rate changes do not produce spurious spectral peaks. We first de-trended the bathymetry and gravity data in the time domain. To help identify spectral peaks for bathymetry, the spectra were calculated on the time-derivative of the bathymetry, which amounts to pre-whitening the data [START_REF] Percival | Spectral Analysis for Physical Applications: Multitaper and Conventional Univariate Techniques[END_REF]. We limited all spectral analysis to wavelengths that occur five or more times in the analyzed section of the geophysical transect (e.g., <5.6

Myr for the two 0-28 Myr regions adjacent to the ridge axis) and were not filtered out via preprocessing (<300 kyr for gravity and <10 kyr for bathymetry). We define a spectral peak as a region that has a higher value than the 95 th percentile confidence interval of the background spectra.

We separated the transect into three sections for analysis: the entire eastern flank (0-28 Ma), the corresponding near ridge western flank (0-28 Ma), and the older, off-axis lithosphere (28-51 Ma) on the western flank. The window of the east flank (blue portion of Figure 6) is set by the length of our transect. The window of the near-axis west flank (red portion of Figure 6) is set equal to the eastern flank for comparison purposes. Further west, the cruise encountered stormy weather that added high frequency noise to our gravity data that the Butterworth filter was unable to remove, so we omit these from the subsequent analysis. Setting these boundaries allows us to employ a standard multi-taper method and intuitively present the analysis and confidence intervals. We find that the two near-axis sections (red and blue, Figure 6) display similar bathymetry and RMBA power spectra. The off-axis section (orange, Figure 6) lacks the same peaks in the bathymetry and has weaker power overall, likely due to mass wasting and sedimentation.

For the North American (western) 0-28 Ma bathymetry power spectra, broad peaks appear at 300-, 390-, 550-, and 950-kyr periods (dashed vertical black lines in Figure 6b). For the African (eastern) bathymetry power spectra, broad peaks appear at 270-, 340-, 550-, and 950-kyr periods. The near-axis RMBA power spectra show peaks at 390-, 550-, and 950-kyr periods on both the North American and African plates. Power is diffuse as spectral densities lack peaks but remain high at longer periods (> 1 Myr) for both the bathymetry and RMBA. The off-axis North American 28-51 Ma RMBA spectral density (orange in Figure 6d) shows a strong peak at 390 kyr and a diffuse peak at 650 kyr (dash-dotted line in Figure 6d).

The shortest period peaks in both the RMBA and bathymetry likely reflect abyssal hill spacing. The typical spacing of abyssal hill-bounding faults along the transect is 2-6 km (Figure 3b). For a half spreading rate of ~14 mm yr -1 , this implies that faulting should add power to the bathymetry spectra between 140-420 kyr. This is consistent with the spectral power in the bathymetry between 250-400 kyr, as well as the 390-kyr spectral peaks in the RMBA, suggesting that faulting modifies the crustal thickness at these periods.

The 550-and 950-kyr peaks in the RMBA of the near-axis sections and the 650kyr diffuse peak in the off-axis section do not correlate with the abyssal hill fabric and may instead reflect processes that result in variations in crustal thickness (see discussion in Section 5). The change in RMBA spectral peaks between the younger, near-axis (0-28 Ma) and the older, off-axis (28-51 Ma) suggests a change in the tectono-magmatic processes causing these peaks. There may also be power at longer wavelengths (e.g., 10-20 Myr) as seen by the periodicity in the raw RMBA data (Figure 3e); however, the section length limits us to interpreting periods in the power spectra less than 5.6 Myr.

There are no significant peaks in the bathymetry power spectrum at any of the Milanković cycle periods (100 kyr, 41 kyr, and 23 kyr) as reported by [START_REF] Crowley | Glacial cycles drive variations in the production of oceanic crust[END_REF], supporting the idea that ocean floor bathymetry formed at slow-spreading MORs is not sensitive to fluctuations in melt supply triggered by sea level change [START_REF] Olive | Sensitivity of seafloor bathymetry to climate-driven fluctuations in mid-ocean ridge magma supply[END_REF]2016;[START_REF] Goff | Comment on "Glacial cycles drive variations in the production of oceanic crust[END_REF][START_REF] Goff | No Evidence for Milankovitch Cycle Influence on Abyssal Hills at Intermediate, Fast, and Superfast Spreading Rates[END_REF].

Admittance & Coherence

Elastic plate thickness, T e, provides an estimate for the thickness of a rigid elastic plate (overlying an inviscid half-space) that best fits the observed patterns of lithospheric flexure. It is therefore a useful proxy for the integrated strength of the plate, and a direct control on the morphology of the seafloor [START_REF] Watts | Isostasy and Flexure of the Lithosphere[END_REF]. To estimate 𝑇 𝑒 , we calculated the admittance and coherence of the topography and MBA along the entire SCARF transect following the approach of Forsyth (1985) (Figure 7). Admittance is defined as the transfer function between topography and MBA. The highest value of admittance is reached through isostatic compensation, which produces topography in response to local thickening of the crust. Short-wavelength oscillations in crustal thickness, however, tend to be flexurally compensated, suppressing strong topographic signals. Thus, admittance is usually low in the short-wavelength limit. Thinner elastic plates will have a larger amplitude admittance at shorter wavelengths because they can isostatically compensate shorter wavelength topography.

Coherence is a statistic that characterizes the strength of the association between two spectra at different wavelengths. For an elastic plate, long wavelength topographic variations will be isostatically compensated, and thus bathymetry and the MBA will have a coherence of 1. On the other hand, short wavelength topography that is not isostatically compensated will have a coherence of 0. For the calculation of admittance and coherence, the bathymetry and MBA data were not pre-whitened and the densities of the crust and mantle were assumed to be 2800 and 3300 kg m -3 , respectively. For comparison, we also calculated theoretical coherence and admittance curves for an idealized elastic plate model of a given 𝑇 𝑒 , assuming a Young's modulus of 100 GPa, Poisson's ratio of 0.25, a crustal thickness of 6 km, and equal loading on the surface and base of the plate. Estimates of 𝑇 𝑒 from coherence are relatively insensitive to changing the assumed ratio of loading on the surface and base of the plate and to reasonable oceanic crustal thickness variations (4-8 km) [START_REF] Forsyth | Subsurface loading and estimates of the flexural rigidity of continental lithosphere[END_REF]. Coherence and admittance estimates of T e are also insensitive to the inclusion or exclusion of the off-axis (28-74 Ma) data.

An elastic plate thickness (𝑇 𝑒 ) of 6-12 km best fits both the admittance and coherence (Figure 7). 𝑇 𝑒 calculated from the admittance and coherence analyses reflects the strength of the plate at the time of loading. Our results are in the range found for 𝑇 𝑒 in the Atlantic based on applying elastic flexure theory to off-axis volcanoes emplaced on young (<10 Ma) oceanic crust [START_REF] Calmant | Anomalous elastic thickness of the oceanic lithosphere in the south-central Pacific[END_REF]. [START_REF] Smith | Fault rotation and core complex formation: Significant processes in seafloor formation at slowspreading mid-ocean ridges (Mid-Atlantic Ridge, 13-15 N)[END_REF] found evidence of highly-rotated low offset faults in the Central Atlantic which would suggest T e < 1 km at ages <1 Myr . Thus, because our admittance and coherence analyses for the entire SCARF transect predicts a 𝑇 𝑒 in agreement with young lithosphere, we interpret that the loading was applied near, but not on the ridge axes, most likely during tectono-magmatic processes associated with the creation of the oceanic lithosphere.

Discussion

Variations in Crustal Thickness

Gravity-derived crustal thickness varies between 3 and 9 km with a standard deviation of 1.0 km assuming a mean of 6 km along the SCARF transect (Figure 3f). These fluctuations are consistent with the range of crustal thicknesses observed in seismic studies in Atlantic crust away from hotspots (e.g., [START_REF] Hooft | Crustal thickness and structure along three contrasting spreading segments of the Mid-Atlantic Ridge, 33.5°-35[END_REF]Canales et al., 2000a,b;Lizzaralde et al., 2004), as well as previous gravity-derived crustal thickness analyses for the Atlantic Ocean [START_REF] Lin | Evidence from gravity data for focused magmatic accretion along the Mid-Atlantic Ridge[END_REF][START_REF] Detrick | Gravity anomalies and crustal thickness variations along the Mid-Atlantic Ridge between 33° N and 40° N[END_REF][START_REF] Wang | Crustal thickness anomalies in the North Atlantic Ocean basin from gravity analysis[END_REF].

We observe spectral peaks in the near-axis RMBA data at 390-, 550-, and 950-kyr periods and diffuse power at periods longer than 1 Myr. These observations are consistent with seismic imaging of oceanic crust, which find crustal thickness variations on timescales of 0.5-1 Myr in the Atlantic (e.g., near the Kane fracture zone: Canales et al., 2000a).

Variations in RMBA-derived crustal thickness could be caused by any one of the following tectonic and magmatic processes: variations at slow to ultra-slow spreading rates, along-segment thickness variations, faults offsetting the Moho, changes in ridge geometry and upper mantle flow, magma solitons, or thermal and chemical variations in the mantle source region (Figure 8). Below we discuss these different processes and their potential influence on the observed crustal thickness variations along the SCARF transect.

Spreading rate is not a dominant cause of crustal thickness variation. Our data show no significant correlation between spreading rate and crustal thickness (Figure A4). This agrees with observations and geodynamic melting models that show crustal thickness is not sensitive to spreading rate for spreading rates greater than 10 mm yr -1 (Parmentier & Morgan, 1990;[START_REF] Chen | Oceanic crustal thickness versus spreading rate[END_REF][START_REF] Behn | Melting systematics in mid-ocean ridge basalts: Application of a plagioclase-spinel melting model to global variations in major element chemistry and crustal thickness[END_REF].

Second, the inferred crustal thickness variations could be due to our transect wandering slightly off the flow line, either closer to the middle or edge of the ridge segment. Indeed, crustal thickness variations between 3 and 9 km have been observed on single ridge segments in the North Atlantic [START_REF] Cannat | How thick is the magmatic crust at slow spreading oceanic ridges?[END_REF][START_REF] Hooft | Crustal thickness and structure along three contrasting spreading segments of the Mid-Atlantic Ridge, 33.5°-35[END_REF]. Thus, if the transect was not always perfectly aligned with the segment center, this could explain some of the variability in our observations. This variation should not have a clear expression in the frequency domain due to its stochastic nature, and thus we do not consider it to be the primary source of crustal thickness variations along the SCARF transect.

Variations in ridge segment geometry and upper mantle flow have also been theorized to affect crustal thickness on both segment and regional scale (e.g., [START_REF] Detrick | Gravity anomalies and crustal thickness variations along the Mid-Atlantic Ridge between 33° N and 40° N[END_REF][START_REF] Hooft | Crustal thickness and structure along three contrasting spreading segments of the Mid-Atlantic Ridge, 33.5°-35[END_REF][START_REF] Bonatti | Mantle thermal pulses below the Mid-Atlantic Ridge and temporal variations in the formation of oceanic lithosphere[END_REF]. Rearrangement and/or changes in the length of segments could therefore bias our one-dimensional transect of crustal thickness observations. The timescales for such changes is typically on the order of millions of years [START_REF] Schouten | The memory of the accreting plate boundary and the continuity of fracture zones[END_REF][START_REF] Caress | Eocene reorganization of the Pacific-Farallon spreading center north of the Mendocino Fracture Zone[END_REF][START_REF] Tucholke | Segmentation and crustal structure of the western Mid-Atlantic Ridge flank, 25° 25′-27° 10′ N and 0-29 my[END_REF]. Segment reorientation is observed through changes in the ridge-parallel magnetic lineation. For example, a lineation-perpendicular offset formed by the growth of a fracture zone or the rotation of the ridge-parallel magnetic lineations due to ridge rotation or rift propagation (e.g. [START_REF] Kleinrock | Fast rift propagation at a slow-spreading ridge[END_REF]. Global magnetic data [START_REF] Meyer | EMAG2: Earth magnetic anomaly grid (2arc-minute resolution) version 3[END_REF] along our trackline show no alteration in ridge-parallel magnetic lineations. Thus, there is no evidence of ridge segment length variability on either side of the ridge axis based on global magnetic data. Furthermore, ridge segment rearrangement is not periodic, so if it did occur, it would not be expected to manifest as discrete spatial peaks.

Normal faulting can locally alter crustal thickness by offsetting the base of the crust. Previous seismic oceanic crustal thickness studies have shown faulting to be a major driver of observed cross-axis crustal thickness variations because crustal thickness anomalies are often found with fault scarps directly above them (e.g., Canales et al., 2000a[START_REF] Seher | Crustal velocity structure of the Lucky Strike segment of the Mid-Atlantic Ridge at 37 N from seismic refraction measurements[END_REF]. As shown in Figure 3b, fault spacing within 400 km of the ridge axis ranges from 2-6 km along the SCARF transect. For an average half spreading rate of ~14 km Myr -1 , the spectral power associated with such faulting should lie in periods of 140-420 kyr, coincident with spectral peaks in the bathymetry data between 250-400 kyr and the RMBA spectrogram peak at 390 kyr. However, the observed faulting periods do not coincide with the additional spectral peaks at 550 and 950 kyr. In the most extreme case, large-offset (>5 km) detachment faults can cause observable crustal thickness variations up to 1-3 km lasting for 1-3 Myr [START_REF] Tucholke | Megamullions and mullion structure defining oceanic metamorphic core complexes on the Mid-Atlantic Ridge[END_REF][START_REF] Parnell-Turner | Genesis of corrugated fault surfaces by strain localization recorded at oceanic detachments[END_REF]. Detachment faults create strong asymmetry in crustal thickness, exposing deep crustal units in the footwall of the fault, while the hanging wall can experience emplacement of a complete crustal section [START_REF] Olive | The structure of oceanic core complexes controlled by the depth distribution of magma emplacement[END_REF]. Because we interpret Kafka Dome as an oceanic core complex in 8 Ma crust east of the ridge axis (see Section 3.6), detachment faulting is a potential factor driving crustal thickness variations. However, the impact of a detachment fault on crustal thickness should be local. In our dataset, this would most likely add diffuse power at periods similar to their lifespan (1-3 Myr) instead of generating clear spectral peaks.

Magma solitons are buoyantly ascending porosity waves of high melt fraction that occur in porous two-phase flow such as the mantle-melting column beneath a mid-ocean ridge [START_REF] Scott | Magma solitons[END_REF]. These solitons could result in periodic crustal thickness variations if they carry sufficient melt to the ridge axis. Compaction and ascent timescales for basaltic melt in the sub-ridge mantle depend on the melting column height, the melt density contrast, the melt viscosity, the melt fraction, the grain size, and the permeability-porosity relationship [START_REF] Mckenzie | The extraction of magma from the crust and mantle[END_REF]. Depending on the choice of parameters, this timescale can range from 10 4 -10 7 years. These timescales reflect the time necessary to extract melt via compaction once, and provide no requirement that the melt extraction process is periodic. More recent two-phase flow models [START_REF] Sim | Melt focusing at mid ocean ridges[END_REF]) at intermediate spreading rates (3.5 cm yr -1 half spreading rate) produce stable and persistent porosity waves that produce changes in crustal thickness with an amplitude of ~0.5 km and recurrence time scales on the order of 100s of kyrs. Further two-phase flow modelling is necessary in order to ascertain whether this periodicity exists at slower spreading rates and over which range of parameter values.

Source region heterogeneities, chemical and thermal, have also been invoked to explain observed variations in oceanic crustal thickness. Mantle compositional heterogeneity can cause variations in crustal production by enriching or depleting the mantle source, which modulates the melting fraction below the ridge (Katz & Weartherley, 2012). Thermal heterogeneities associated with mantle hotspots change the melt flux of the upwelling mantle. [START_REF] Bonatti | Mantle thermal pulses below the Mid-Atlantic Ridge and temporal variations in the formation of oceanic lithosphere[END_REF] interpreted a combination of gravity and geochemical data collected along the Vema fracture zone as indicating an increase in mantle potential temperature in this region of the Central Atlantic Ocean over the last 20 Myr, causing an increase in crustal thickness on the timescale of the temperature variations (>10 Myr). However, regional gravity-derived crustal thickness calculations [START_REF] Wang | Crustal thickness anomalies in the North Atlantic Ocean basin from gravity analysis[END_REF] suggest more complicated trends in this region.

Plume thermal pulsing has been interpreted to create regional crustal thickness variations that are reflected in the V-shaped ridges near the Azores. However, these Vshaped ridges do not extend south to the S. AMAR segment and [START_REF] Escartín | Crustal thickness of V-shaped ridges south of the Azores: Interaction of the Mid-Atlantic Ridge (36-39 N) and the Azores hot spot[END_REF] argue that the melt flux anomaly from the Azores plume stops at approximately 36˚ N, ~30 km north of the SCARF transect. Thus, while thermal pulsing can affect crustal thickness, there is no strong evidence for thermal source heterogeneities along the SCARF transect.

Trace element and isotopic data in dredged gabbros suggest chemical heterogeneities near the SCARF transect caused by influence from the Azores plume or patches of former sub-continental mantle lithosphere [START_REF] Bougault | Mantle heterogeneity from trace elements: MAR triple junction near 14 N[END_REF][START_REF] Shirey | Three-component isotopic heterogeneity near the Oceanographer transform, Mid-Atlantic Ridge[END_REF]. To cause temporal variations of crustal thickness, these source heterogeneities would have to pass through the melting column entirely; thus the timescales of variation would depend on their size and spacing given a constant mantle upwelling velocity. For instance, in order for the 550-and 950-kyr spectral peaks to correspond to melting anomalies associated with mantle source heterogeneities upwelling at a rate of ~14 km Myr -1 , these heterogeneities would need to have a characteristic size and/or spacing of 7.7 and 13.3 km. Larger heterogeneities would cause longer period crustal thickness variations.

Because we do not have a long time-series chemical dataset near this ridge-segment, we lack the ability to determine how much of the crustal thickness variations are due to chemical source heterogeneities.

In summary, diffuse power at longer periods (>1 Myr) along the SCARF transect is most likely due to a combination of detachment faulting, mantle source heterogeneities (thermal and/or compositional), and/or magma solitons. The lack of peaks at long periods is likely due to the lack of periodicity and the temporal independence of these processes, allowing the expression of their signals in crustal thickness to constructively or destructively interfere. Faulting is most likely the cause for spectral power in bathymetry and RMBA between 250-400 kyr because this time scale matches the characteristic spacing of fault-bounded abyssal hills. By contrast, only regularly spaced chemical source heterogeneities and/or magma solitons can explain the spectral peaks at 550 and 950 kyr. Further numerical studies are necessary to help estimate the relative importance and plausibility of the different factors contributing to the 550-and 950-kyr periods. The change in RMBA spectral peaks from 550 and 950 kyr in the younger, near-axis sections (0-28 Ma) to a diffuse 650-kyr peak in the older, off-axis section (28-51 Ma) could be caused by a change in these tectono-magmatic processes.

Variations of Faulting Style and Relation to Magma Supply Variations

Geophysical modeling shows that the main control on fault spacing is the fraction of spreading accommodated by magmatic accretion, M [START_REF] Buck | Modes of faulting at mid-ocean ridges[END_REF]Behn & Ito, 2008;Ito & Behn, 2008;[START_REF] Tucholke | Role of melt supply in oceanic detachment faulting and formation of megamullions[END_REF][START_REF] Olive | Sensitivity of seafloor bathymetry to climate-driven fluctuations in mid-ocean ridge magma supply[END_REF][START_REF] Liu | Magmatic controls on axial relief and faulting at midocean ridges[END_REF][START_REF] Howell | Seafloor expression of oceanic detachment faulting reflects gradients in mid-ocean ridge magma supply[END_REF]. Co-located analyses of ocean floor bathymetry and basalt geochemistry have validated the importance of magma supply on abyssal hill morphology [START_REF] Roth | Discrete characterization of abyssal hills: Unraveling temporal variations in crustal accretion processes at the 10° 30′ N segment, East Pacific Rise[END_REF]. With higher M, faults are advected off axis into thicker lithosphere more rapidly. This makes continued slip on the existing fault energetically less favorable [START_REF] Forsyth | Finite extension and low-angle normal faulting[END_REF][START_REF] Buck | Effect of lithospheric thickness on the formation of high-and lowangle normal faults[END_REF], resulting in the formation of a new fault near the ridge axis, and an overall decrease in fault spacing and fault offset (Behn & Ito, 2008). A special case occurs when M approaches 0.5, such that spreading on one side of the ridge is accommodated purely by magma intrusion, while extension on the conjugate side is entirely accommodated by slip on a long-lived detachment fault whose footwall forms an oceanic core complex [START_REF] Buck | Modes of faulting at mid-ocean ridges[END_REF]. These models predict an inverse relationship between fault spacing and crustal thickness.

These model predictions imply that if Kafka Dome is indeed the footwall of a detachment fault, the minimum M observed should be close to 0.5 (e.g., [START_REF] Macleod | Life cycle of oceanic core complexes[END_REF][START_REF] Schouten | Tectonic versus magmatic extension in the presence of core complexes at slow-spreading ridges from a visualization of faulted seafloor topography[END_REF]. The maximum M calculated through fault heave analysis is 0.85 (Figure 3c). This range is in agreement with the overall range in M values (0.45-0.75) that has been calculated in young lithosphere along the MAR (Behn & Ito, 2008[START_REF] Paulatto | Heterogeneous and asymmetric crustal accretion: New constraints from multibeam bathymetry and potential field data from the Rainbow area of the Mid-Atlantic Ridge (36˚ 15'N)[END_REF].

The geophysical data collected along the SCARF transect allows us to further quantify the relationship between magma supply and fault spacing. The Spearman rank correlation coefficient [START_REF] Spearman | The proof and measurement of association between two things[END_REF], r, detects any type of monotonic correlation rather than only a specific functional correlation and is less sensitive to outliers. Values of r range from -1 to 1, with larger absolute values indicating that the two variables more strongly co-vary according to a monotonically increasing (positive) or decreasing (negative) relationship. The corresponding p-value indicates the probability that the relationship is due to randomness. We find a strong and significant negative correlation (r =-0.64, p=2.2e-4) between average crustal thickness and fault spacing within 400 km of the ridge (Figure 9a). This agrees with previous findings in the Atlantic Ocean of fault spacing increasing at the edge of ridge segments where the crust is thinner [START_REF] Shaw | Ridge segmentation, faulting and crustal thickness in the Atlantic Ocean[END_REF][START_REF] Shaw | Causes and consequences of variations in faulting style at the Mid-Atlantic Ridge[END_REF]. Thus, we interpret the variations in fault spacing along the flow line (Figure 3b) to reflect variability in magmatic input into the brittle crust of the MAR over time. We find no significant correlation between the calculated M value and crustal thickness over the same bins (r=-0.254, p=0.19) (Figure 9b). This lack of a significant correlation could be due to the coarse resolution of our calculated M values, inaccuracies in our heave measurements due to off-axis sedimentation and mass wasting, or that M is independent of crustal thickness as observed by [START_REF] Howell | Magmatic and tectonic extension at the Chile Ridge: Evidence for mantle controls on ridge segmentation[END_REF]. [START_REF] Olive | The structure of oceanic core complexes controlled by the depth distribution of magma emplacement[END_REF] modeled extensional faulting at mid-ocean ridges with separated melt accretion rates in the ductile and brittle layers of oceanic crust and showed that faulting is only sensitive to the M value in the brittle lithosphere. Thus, varying amounts of magma intrusion in the asthenosphere does not affect faulting style and could result in a decoupling of fault style and crustal thickness variations.

Finally, we note that r statistics can be sensitive to outliers. For example, removing the two data points with high M and high crustal thickness results in a nonintuitive negative correlation (r=-0.495, p=0.01), illustrating the weakness of the correlation between M and crustal thickness. Thus, while a correlation between M and crustal thickness is expected on a global scale-as both are known to increase with increasing spreading rate [e.g., Dick et al., 2003;[START_REF] Olive | Sensitivity of seafloor bathymetry to climate-driven fluctuations in mid-ocean ridge magma supply[END_REF]-it is possible that our transect did not encounter strong enough fluctuations in magma supply to be clearly reflected in both crustal thickness and tectonic fabric.

Conclusion:

We conducted one of the longest continuous geophysical surveys along a 1400km spreading-parallel flow line across the Mid-Atlantic Ridge from 28 Ma on the African Plate to 74 Ma on the North American plate. The transect was analyzed to study covariability in bathymetry, gravity, and fault statistics in order to elucidate tectonomagmatic variability through time at a single slow-spreading ridge segment. Gravityderived crustal thickness varies from 3-9 km along the transect with a standard deviation of 1.0 km assuming a mean crustal thickness of 6 km. Admittance and coherence of gravity and bathymetry predict the effective elastic thickness of the lithosphere, 𝑇 𝑒 , to be 6-12 km, consistent with other estimates for young oceanic lithosphere and indicating that most lithospheric loading that produced relief occurred near, but not on the ridge axis.

Bathymetry and RMBA spectral densities show concurrent peaks on both the eastern and western flanks (0-28 Ma) at 390-, 550-, and 950-kyr periods, with diffuse power at longer periods (>1 Myr). We interpret the diffuse power along the transect as due to a combination of detachment faulting, mantle source heterogeneities (thermal and/or compositional), and variations in upper mantle flow. The 390-kyr spectral peak corresponds to the characteristic spacing of abyssal hill-bounding faults. By contrast, the 550-and 950-kyr peaks are more likely explained by short-wavelength mantle source heterogeneities and/or magma solitons. Further off-axis (28-51 Ma), RMBA spectral density shows a diffuse peak at 650 kyr instead of 550 and 950 kyr could be caused by a change in these tectono-magmatic processes over >10 Myr time periods. Further modelling and analysis are necessary to better quantify the characteristic time scales of these processes and their parameter dependence.

Finally, fault statistics and the identification of a newly observed, off-axis oceanic core complex suggest M values ranging from 0.5-0.85 along flow line. We find a significant negative correlation between fault spacing and gravity-derived crustal thickness, indicating the variations in fault spacing along the flow line reflect variability in magmatic input into the brittle crust of the MAR over time. This implies that future statistical analyses of fault populations have the potential to identify variations in mantle melting and melt transport at mid-ocean ridges. In particular, further studies that jointly analyze faulting characteristics and crustal thickness have the potential to elucidate past variations in magma production.

Acknowledgements:

Data and supplemental materials are available at the Woods Hole Open Access Server (doi.org/10.26025/1912/24796). We would like to thank the Woods Hole Oceanographic Institution, National Science Foundation, Naval Oceanographic Office, and the captain and crew of R/V Neil Armstrong for making the SCARF cruise possible.

We would also like to thank Eboné Pierce for her help during the cruise. We thank Meghan Jones for advice using MBSystem. We also thank Maurice Tivey, John Greene, Figure Captions: Figure 1: Comparison of shipboard and regional data sets along the SCARF cruise track. a) Track line bathymetry (scale in km is shown by white bar) plotted along the track line (straight black line) over satellite-derived bathymetry [START_REF] Weatherall | A new digital bathymetric model of the world's oceans[END_REF] for the region. The inset shows global plate boundaries (black lines) and our track line (red line). b) Shipboard magnetic anomaly (scale in nT is shown by white bar) plotted along the track line over global magnetic anomaly grid EMAG2v3 [START_REF] Meyer | EMAG2: Earth magnetic anomaly grid (2arc-minute resolution) version 3[END_REF]. White lines are other shipboard magnetic anomaly tracks from NCEI's Marine Trackline Geophysical database (https://www.ngdc.noaa.gov/mgg/geodas/trackline.html). c) Shipboard free air gravity anomaly (scale in mGal is shown by white bar) plotted along the track line against 1 arcminute satellite-derived global free air anomaly version 24 [START_REF] Sandwell | New global marine gravity model from CryoSat-2 and Jason-1 reveals buried tectonic structure[END_REF]. Figure 2: Comparison of the observed magnetic anomaly along the SCARF transect (blue line) and synthetic magnetic anomaly (black line) calculated from the magnetostratigraphic scale plotted on the bottom. See section 3.1 for synthetic calculation method. Synthetic profile is vertically translated by -700 nT for illustration purposes. Isochron picks are labeled and represented as dashed grey lines. Figure 3: Shipboard geophysical datasets and derived geophysical parameters along the SCARF transect. a) Track line bathymetry (black) with picked fault locations (red dots, displayed 1 km shallower than picked location). The location of the observed oceanic core complex, Kafka Dome is marked by the blue star, displayed 0.7 km deeper than the picked location. b) Binned fault spacing (open circles) and interpolated sediment thickness from compiled isopach maps, ocean drilling records, and seismic reflection data (red line; [START_REF] Divins | Total sediment thickness of the world's oceans & marginal seas[END_REF]. Grey regions denote potential oceanic core complex regimes located symmetrically across the ridge axis. c) M (circles) calculated as described in Section 3.2. d) Half spreading rates calculated from our magnetic anomaly derived ages (black line) and from an interpolated global age grid [START_REF] Müller | Age, spreading rates, and spreading asymmetry of the world's ocean crust[END_REF]. e) Calculated residual mantle Bouguer anomaly (RMBA). f) Downward continued crustal thickness variations calculated from the RMBA assuming an average crustal thickness of 6 km. Figure 4: Cumulative frequency distribution of fault throw for the entire SCARF transect (solid black line), as well as for only those faults within 200 km (blue) and 400 km (red) of the ridge axis. The throw distribution from near-ridge picked faults on the Chile Ridge (dashed black line) is shown for comparison. Straight lines in this plot indicate an exponential fault distribution, which is predicted for natural normal fault populations (see Section 3.2). We interpret the difference in the fault population between those located within 400 km of the ridge to that calculated for the entire transect to be due to the increasing effects of sedimentation and mass wasting. Thus, we limit our interpretations to fault statistics within 400 km. The change in slope of the fault distributions for throws > 0.6 km is likely due to the low number of observed high throw faults in our dataset. Figure 5: a) Shaded bathymetry for our observed oceanic core complex, Kafka Dome. White line denotes the ship track. b) Calculated RMBA along the track line showing an 8 mGal local increase over the Kafka Dome, characteristic of many oceanic core complexes. c) Bathymetry from oceanic core complex 13˚20' N [START_REF] Escartín | Tectonic structure, evolution, and the nature of oceanic core complexes and their detachment fault zones (13 20' N and 13 30' N, Mid Atlantic Ridge)[END_REF] for comparison. The horizontal scale is the same in both a) and c).

Figure 6 :

 6 a) Track line bathymetry separated into four sections. b) Power spectral density of pre-whitened bathymetry for the three colored sections. c) Calculated residual mantle Bouguer anomaly (RMBA) separated into the same sections. d) Power spectral density of pre-whitened RMBA. The approximate 5 th to 95 th percentile confidence interval from the spectral method is indicated by the bottom and top of the black bar in the upper right corner of b) and d) relative to a given value (circle). Vertical dashed lines show interpreted concurrent peaks between the east (red) and the near-axis west (blue) at 390-, 550-, and 950-kyr periods. Dash-dotted line in d) shows the interpreted diffuse peak at 650 kyr for the 28-51 Ma west flank RMBA. Off-axis data shows a different RMBA and bathymetry spectrum compared to the near ridge data.

Figure 7 :

 7 a) Calculated admittance (black line) plotted alongside theoretical admittance curves for different values of effective elastic thickness, 𝑇 𝑒 . b) Calculated coherence (black line) plotted alongside theoretical coherence curves for various values of 𝑇 𝑒 . Both models predict 𝑇 𝑒 values between 6-12 km.

Figure 8 :

 8 Schematic diagram depicting potential causes for oceanic crustal thickness variation discussed in Section 5.1.

Figure 9 :

 9 a) Binned fault spacing versus binned crustal thickness for all bins less than 400 km from the ridge. b) Binned crustal thickness against M for all bins within 400 km of the ridge. Overlapping bins are 100 km wide every 25 km.

  Figure 1.

and Masako Tominaga for advice on processing the magnetic datasets. We would like to thank Peter Huybers for sharing his spectral analysis codes. We would like to thank Rob Sohn for his help on interpreting spectral analysis. We would like to thank Del Bohnenstiel, Milena Marjanović, one anonymous reviewer, and Editor Thorsten Becker for their very helpful comments that improved this manuscript. Funding was provided for this research by NSF OCE-14-58201.


Francesco Scotognella
email: francesco.scotognella@polimi.it

Giuseppe Della Valle

Ajay Ram

Srimath Kandada

Dirk Dorfs

Margherita Zavelani-Rossi

Matteo Conforti

Guglielmo Lanzani

Liberato Manna

Francesco Tassone

Plasmon dynamics in colloidal Cu 2-x Se nanocrystals

Theoretical and experimental research on the optical features exhibited by metallic nano-systems has paved the way to extensive applications in several fields, from surface-enhanced spectroscopies, to biological and chemical nano-sensing. Such developments are primarily due to the localized plasmon resonances (LPR), resulting in intense optical absorption and scattering as well as sub-wavelength localization of large electrical fields in the vicinity of the nano-object. [START_REF] Schmid | Clusters and colloids: from theory to applications[END_REF][START_REF] Tam | [END_REF][3][START_REF]Metal nanoparticles: synthesis, characterization, and applications[END_REF][START_REF] Zeman | [END_REF][6][7][8][9] The plasmonic response of these metallic nanostructures is strongly dependent on the type of metal of which they are made, on the dielectric function of the surrounding medium, on the particle shape and, to a lesser extent, it is also dependent on the particle size. However, for ultra-small nanoparticles, size dependent effects are negligible, and for homogeneous (i.e. purely metallic) structures, the plasmonic resonance is only determined by the properties of the metal. Since few metallic materials such as Ag, Au, Cu and Pt have been successfully shaped into ultra-small nanoparticles, [10][11][12][13][14] plasmonic response has been mostly reported in the visible, with the exception of some Au-based elongated nanostructures which have also demonstrated plasmonic response in the near infra-red (~800 nm). [15][16] In order to circumvent such a limitation, a consistent effort has been spent in the last years in the synthesis of heavily-doped semiconductor nanocrystals so as to achieve metallic behaviour and eventually a tunable plasmonic response in a broader optical spectrum. In particular, research was concentrated onto S-based, Se-based and Te-based II-VI semiconductor nanocrystals with Cd-Cu substitution. [17][18][19][20][21][22] In these nanoparticles Copper can exist in a wide range of stoichiometric ratios in addition to the 2:1 with respect to Selenium or Sulfur, resulting into Copper vacancies and into "selfdoping" of the material. The influence of this heavy doping on the optical response in the infrared region has been a subject of interest in the last decade. In 1973, Gorbachev and Putilin observed a plasmon band in the reflectivity of p-type copper selenide and copper telluride thin films. [23] They also found that the spectral position of this band strongly depends on temperature, a fact that they attributed to a strong dependence of the hole effective mass on temperature. Very recently, Alivisatos and coworkers [22] and Manna and coworkers [24] reported a detailed study on plasmonic properties of Cu 2 X nanocrystals (with S in [22] and Se in [24], respectively), which showed how starting from nearly stoichiometric Cu 2 X particles, their oxidation (upon exposure on air) increased the degree of Copper deficiency leading to the emergence of an IR band that gradually increased in intensity and blue-shifted in energy. Manna et al. also showed the reversible tuning of the IR plasmonic response of Cu 2-x Se nanocrystals via their gradual oxidation, either under air or by stepwise addition of a Ce(IV) complex. [24] In the very last years, the dynamical response of plasmonic nano-structures has also attracted more and more attention, giving rise to the emerging field of ultra-fast active plasmonics [25]. Actually, it has been recently reported that direct optical excitation of the metal by intense femtosecond laser pulses lead to ultra-fast modulation of plasmonic signals, thus resulting in unprecedented terahertz modulation bandwidth, a speed at least five orders of magnitude faster than existing technologies. [25] A quantitative investigation of such dynamical features is thus of the most importance for the development of next generation of ultra-fast nano-devices. Pump-probe spectroscopy, giving access to the electronic excitation and subsequent relaxation processes in the material, has been demonstrated to be the eligible tool for the experimental study of the ultra-fast dynamical features exhibited by plasmonic nanostructures. So far, several noble metal sctructures have been investigated, including spherical particles [26][27][28], as well as nanorods [29], leading to a general scenario of the underlying physical phenomena similar to the one observed in bulk (i.e. thin film) metallic systems [30,31]. Actually, when the pump pulse is in the infrared the initial Fermi distribution of electrons in the conduction band is strongly perturbed by pump absorption; the pump pulse creates energetic electrons that are not in thermal equilibrium and within few hundred fs a new Fermi distribution is achieved by strong electronelectron scattering, resulting in a thermalized electron gas with higher temperature than the lattice (hot electrons); subsequently, within the next few ps the electron gas cools down by releasing its excess energy to the lattice through electron-phonon coupling; and, ultimately, within hundred of ps, the nanoparticle releases its energy to the environment, [16] (with heat conduction to the surface of the nanoparticles provided by phonon-phonon coupling).

In this work, we report on the optical response of Cu 2-x Se colloidal nanocrystals to ultrafast laser pulses, recorded via pump-probe experiments in the infra-red with 200 fs resolution time. In particular we investigated the hole-phonon coupling dynamics taking place in the ps time scale. The experimental results are interpreted according to the two-temperature model already developed for the theoretical investigation of more conventional metallic systems.

The Cu 2-x Se nanocrystals considered in the present study have been prepared following the synthesis process reported in [24]. As already demonstrated, the optical response exhibited by these nanostructures under continuous wave excitation with infra-red light is plasmonic in nature and can be quantitatively interpreted according to quasi-static approximation of the Mie theory [24]. More precisely, the dipolar absorption and scattering cross-sections of a Cu 2-x Se nanoparticle of radius R << 2πc/ω are respectively given by: with ε ∞ the high frequency dielectric permittivity (due to the presence of inter-band transitions at a higher photon energy), ω P the plasma frequency of the free carriers of the system, and Γ the free carrier damping (i.e. the inverse of the carrier relaxation time).

Our previous results [24] showed that, depending on the 'oxidation' parameter x, different material dielectric functions ε are retrieved, resulting in a plasmon resonance (at frequency ω 0 approximately given by ε 1 (ω 0 ) = -2ε m ) that blue-shifts with increasing x. For the present study we selected x = 0.15 (i.e. Cu 1.85 Se) samples leading to a large blue shift of the plasmon peak. We measured a relatively narrow plasmonic resonance, with an intense peak around 1050 nm, is measured for Cu 1.85 Se particles dispersed in toluene solution (refractive index n m = 1.505, ε m = 2.25), as illustrated in Figure 1, where we plot the optical density measured in a 1 mm thick quartz cuvette with an estimated particle concentration N = 3x10 13 cm -3 . Figure 1 also reports the total extinction cross-section σ E = σ A +σ S of Cu 1.85 Se nanoparticles computed from Eqs. ( 1)-( 2) with Drude parameters ε ∞ = 1.38, ω P = 4.43x10 15 rad/s, and Γ = 9.72x10 14 rad/s, in agreement with our previous study [24].

The pump probe set-up is based on a commercial Ti:Sapphire amplified laser system delivering 150 fs pulses at 1 KHz repetition rate at a central wavelength of 800 nm. A fraction of the beam is used to pump a non-collinear parametric amplifier (NOPA) to generate pulses in the infra-red [32,33]; we tuned the NOPA to obtain pump pulses in correspondence with the plasmonic resonance of the sample (∼1050 nm, 1.19 eV with a bandwidth of 40 nm). The probe pulses are produced by focusing the fundamental beam into a 2 mm thick sapphire plate to generate a stable white light super-continuum. A long pass filter with cut-on wavelength at 820 nm is used to filter out the residual fundamental and the visible component of the probe pulses. The pump and probe beams are focused onto the sample with a spot size of 200 µm. The pump-probe setup, employs a computer controlled optical multichannel analyser and the measured signal is a map of the chirp-free differential transmission ΔT/T = (T on -T off ) / T off as a function of wavelength and the pump-probe time delay; T on and T off are the probe spectra transmitted by the excited and unperturbed sample, respectively. Figure 2 reports the temporal dynamics of the relative differential transmission (ΔT/T) probed at 900 nm, obtained by exciting the nanocrystals with different pump fluencies. The experimental results reveal ΔT/T > 0 right after the absorption of the pump beam, with maximum value of about 40% under the maximum pump fluence of 2.23 mJ/cm 2 . A monothonic and fast decrease of the signal is then observed within a few ps (Figure 2a) followed by a much slower (ns time scale) decay (Figure 2b), leading to complete recovery of the initial condition (before pump arrival) within few ns.

To tentatively describe the non-linear optical features as reported in Fig. 2 we applied the theoretical methods developed for metals, since the linear optical response of Cu 2-x Se nanocrystals is dominated by the metallic behaviour of the system, which leads to a plasmonic resonance in the extinction spectrum, we applied the theoretical methods developed for metals to tentatively describe the non-linear optical features as reported in Figure 2.

The pump-probe dynamics of metallic systems can be ascribed to the variations attained by the metal dielectric function induced by pump absorption. [30,31] As explained before, the ps dynamics is related to the cooling of the hot electron gas generated by absorption of the pump energy. Therefore, the time evolution of the system wihin ps time scale can be modelled as a heat transfer between a thermalized gas of carriers at temperature T C and the lattice at a lower temperature T L , via carrier-phonon scattering process. The so-called Two Temperature Model (TTM) quantitatively accounts for such transfer according to the following coupled equations [30]:

(3)
where γT C is the heat capacity of the carrier gas, being γ the so-called carrier heat capacity constant, C L is the heat capacity of the lattice, G is the carrier-lattice coupling factor, G L is the lattice-environment coupling factor, T 0 is the environmental temperature (which is supposed to be constant) and P A (t) is the pump power density absorbed in the volume of the metallic system. Once the temperature dynamics induced by the pump pulse is determined from Eq. ( 3), the optical response of the systems to the probe pulse can then be retrieved from the temperature dependence of the dielectric function of the material and from the electromagnetic boundary conditions due to the particular geometry of the system, which in the case of a small spherical particle are described by Eq. ( 1).

The effects of the carrier and lattice temperatures on the dielectric function also depend on the detailed band structure of the material, and on the energy of the probe photon. Little is known on the band structure of Cu 2-x Se, apart from the fact that it behaves as a p-type degenerate semiconductor with half filled valence band. [34,24] Moreover, many possible crystallographic phases of Cu 2-x Se have also been found. [35] Here we report in Figure 3.a a tentative and qualitative sketch of the band structure deduced from spectroscopical investigations reported in the literature. According to Al-Mamun et. al. [36], the Cu 2-x Se system with x = 0.15 exhibits two interband optical transitions, a direct transition at 2.1 eV and an indirect transition at 1.3 eV. It is thus expected that the pump photon at 1.19 eV energy is absorbed by an intra-band process similarly to what happens in noble metals in the infra-red. The subsequent heating of the carrier gas results in a smearing of the Fermi distribution, and gives rise to a modulation of the inter-band transition probability for the probe light at 1.38 eV (Figure 3.a), which results in a modulation Δε ∞ of the ε ∞ parameter in the Drude dielectric function Eq. ( 2). Being the probe photon energy here used much lower than the direct inter-band transition edge, and considering that the indirect inter-band process, and its modulation, is much weaker compared to direct transitions, we conclude that the imaginary part of Δε ∞ is negligible compared to the real part. In metals, it is well known that this modulation is proportional to the carrier excess energy [28,31], and thus it scales quadratically with the carrier temperature T C , i.e. Δε ∞ = δ T C 2 , δ being a fitting parameter.

Not only the carrier temperature but also the lattice temperature results in a modulation of the dielectric function ε ∞ . In metallic systems the Drude damping Γ factor increases linearly with T L . However, the resulting temperature dependence of the dielectric function ε is relatively weak and tipically negligible, unless the lattice temperature is increased significantly. Instead Cu 2-x Se exhibits the additional and distinctive dependence of the hole mass on the lattice temperature as reported by Gorbachev and Putilin.

[23] This dependence on the temperature is much stronger compared to the modulation of Γ. For the aims of the present work, we retrieved an analytical expression for such dependence by simply fitting the data reported in [23], providing the following formula for the hole effective mass: m h (T L ) = (aT L b + c) m 0 , with a = 4.04x10 -16 (K -b ), b = 5.71, c = 0.151, and m 0 the free electron mass. The temperature dependence of m h results in a temperature dependence of the plasma frequency of the form:

ω P (T L ) = [N C e 2 /m h (T L )ε 0 ] 1/2 (4)
Therefore the plasma frequency decreases as the lattice is heated (N C being the free carrier density). The effect on the extinction cross-section, numerically computed from quasi-static formulas given in Eq. ( 1) is illustrated in Figure 3b. Note that for the probe wavelength of 900 nm the extinction is strongly reduced as the lattice is heated by just few tens of degrees.

The previous discussion provides a qualitative understanding of the dynamics observed in the experimental data reported in Figure 2.a: the T C dependent inter-band effect, resulting in an increase of the real part of the dielectric function ε of the material, leads to a red-shift of the plasmonic resonance (which is now provided by ε 1 (ω 0 ) + Δε ∞ = -2ε m ), and thus enhances the transmissivity of the probe (which is in the blue-tail of the extinction spectrum, Figure 1) right after the pump absorpion; the T L dependent intra-band effect, due to the temperature dependence of the free carrier effective mass, also results in a red shift of the plasmonic resonance, and is related to the long-living differential transmission tail which can be observed after the end of the carrier-lattice interaction, and for several ps.

The long time-scale dynamics (Figure 2.b) is thus completely decoupled from the fast carrier-phonon dynamics, and thus an exponential decay with a time constant 1/G L is expected.

In order to attempt a quantitative comparison with the experiment, the theoretical differential transmission was computed as ΔT/T = exp(-Δσ E NL) -1, Δσ E being the variation attained by the extinction cross-section according to the quasi-static formulas of Eq.( 1), and the temperature-dependent dielectric function of Cu 2-x Se given by:

(

with T C and T L provided by numerical solution of the TTM, and G, δ, and G L for the long time-scale dynamics, as fitting parameters. In doing so, we assumed a heat capacity of the lattice C L = 2.72x10 6 J m -3 K -1 , according to the experimental data reported in [35] for α-Cu 2 Se, whereas the heat capacity of the carriers was estimated as γ = rγ Au = 1.37 J m -3 K -2 , with γ Au = 62.9 J m -3 K -2 being the heat capacity of the free electrons in Gold and r the ratio between Cu 2 Se carrier density (estimated, from the plasma frequency and hole effective mass at room temperature, to be N C = ω P 2 ε 0 m h (T C )/e 2 = 1.28x10 27 cm -3 ) and gold carrier density (5.9x10 28 cm -3 ).

Given the incident pump fluences of 96, 318 and 2230 µJ/cm 2 , the absorbed pump power density P A (t)

was estimated from pump absorption measurements, showing that only 1% of the incident pump energy is absorbed by the sample, despite the pump wavelength is at the peak of the linear absorption of Cu 1.85 Se nanoparticles. This indicates that the high fluence of the pump leads to strongly non-linear absorption phenomena even during the pump pulse. Such phenomena might be related to the strong perturbation of the plasmonic regime due to energetic and non-thermalized carriers, which could lead to a substantial shift of the plasmonic resonance, and/or to a substantial smearing, so as to drastically reduce the absorption cross-section at the pump wavelength. The absorption dynamics is thus expected to be dominated by the non-thermalized carrier dynamics, which occurs in the tens of fs time scale and is far beyond the resolution time of the present experiment.

Numerical solution of the TTM model with the best fitting parameters is reported in Figure 4a, whereas the corresponding numerically computed differential transmission is shown in Figure 2.a (for the ps time scale) and in Figure 2.b (for the long-time scale). We estimated a low temperature carrier-phonon coupling factor G = G 0 = 2.0x10 15 W m -3 K -1 , which is one order of magniture lower than in Gold. Since the G 0 factor is proportional to the carrier density [37], this lower G is in agreement with the much lower carrier density in Cu 1.85 Se compared to gold. We also estimated δ = 5x10 -8 K -2 , which is about 40% higher than in silver structures at 900 nm probe wavelength [31].

Note that numerical simulations are in good agreement with the experimental data for the lower pump fluences whereas the higher fluence results in some quantitative misagreement within the short timescale dynamics (compare solid and dotted green curves in Figure 2a). Actually for such a high fluence the TTM predicts a carrier temperature increase exceeding 3000 K (see green curves in Figure 4a). It is thus expected that, similarly to what happens in noble metals, the carrier heat capacity γ and the electron-phonon coupling factor G are no more temperature independent [37]. We found that by tentatively introducing temperature dependences for γ and G that mimick the temperature dependences reported for Gold at high temperatures [37], i.e. γ =γ 0 +β 1 (T C -T G ) and G = G 0 +β 2 (T C -T G ) 3 for T C > T G , with β 1 , β 2 and T G fitting parameters, a much better agreement between the model and the experiment is found also for the higher pump fluence (see dashed green curve in Figure 2.a, and see Figure 4.c for fitted temperature dependent parameters). However, being γ and G strongly dependent on the actual DOS of the free carriers [37], a quantitative determination of their temperature dependence can not be asserted without a better knowledge of the band structure of the Cu 1.85 Se crystal, which is out of the scope of the present work.

In conclusion, we report the first observation of the ultrafast optical response of Cu 2-x Se nanoparticles excited at the plasmon resonant energy. We have been able to quantitatively describe the optical response, and to extract from it the electron-phonon coupling constant, which compares reasonably well with those found in noble metals after accounting for the lower carrier density. We remark that the relatively low carrier density and structural peculiarities of Cu 2-x Se lead to significant deviations of the response from the well-known behaviour found in the more traditional metallic nanoparticles (for example the noble metal nanoparticles). In particular, the lower carrier density is responsible for a lower carrier heat capacity and for a much larger effective carrier temperatures at comparable fluences. This leads to strong non-linear effects at high fluences both during the pump absorption process, and during the carrier temperature relaxation process. Moreover, the large temperature dependence of the hole effective mass results into large residual differential transmission signals due to significant shifts of the plasmonic resonance. Such a marked temperature dependence of the hole mass is possibly related to the complex phase diagram of Cu 2-x Se [35], to incipient changes of the crystalline phase, and to a strong sensitivity of the valence band shape on this change of phase. From a technological point of view, this material can be promising, as an example, for the development of ultrafast nano-devices operating at terahertz modulation bandwidth, [25] with the peculiar advantage of plasmon resonance tunability in the IR region, with potential interest for telecom applications.

( 1 )

 1 where ε m = n m 2 is the dielectric constant of the environmental medium, k = n m ω/c and ε = ε 1 (ω) + iε 2 (ω) is the Cu 2-x Se dielectric function provided by the Drude model: (2)

Figure 1 .

 1 Figure 1. Experimental optical density spectrum of the Cu 1.85 Se nanoparticles dispersed in toluene

Figure 2 .

 2 Figure 2. Differential transmission signal from Cu 2-x Se nanocrystal dispersion in toluene for (a) short-

Figure 3 .

 3 Figure 3. (a) Tentative sketch of Cu 2-x Se band diagram for optical transitions in the visible. Smearing of

Figure 4 . 3 K - 1 ;

 431 Figure 4. (a) Carrier (hole) and lattice temperature numerically computed from the TTM of Eq. (3). (b)