High-Resolution Fault Geometry and Motion of the August 23, 2011, M_w 5.8 Mineral, Virginia, Aftershock Sequence Revealed Using Double-Difference and Moment-Tensor Inversion Methods Open Access

The M_W 5.8 Mineral, Virginia, earthquake of August 23, 2011, is one of the largest seismic events known to have occurred east of the Rocky Mountains in historic times and the largest known earthquake in the Central Virginia Seismic Zone (CVSZ). The earthquake caused significant structural damage in Washington, D.C. and prompted a temporary shutdown of the North Anna Nuclear Power Plant, located only 22 km from the epicenter (Figure 1). Many mapped faults exist in the CVSZ, but the Mineral, Virginia, earthquake has not been associated with a known fault structure [1, 2] (Figure 1). Previous studies have determined the regional crustal structure, fault locations, and geologic history of the CVSZ [3-15] while others have attempted to attribute the Mineral, Virginia earthquake to a known fault structure or a previously unmapped fault structure [2, 16-27]. The objectives of this study were to use arrival-time data and P- and S-waveforms of the M_W 5.8 Mineral, Virginia, earthquake aftershock sequence to calculate precise aftershock locations in order to ascertain a high-resolution model of the fault structure illuminated by the aftershock sequence and to evaluate focal mechanisms for aftershocks at various locations on the fault structure to determine fault structure geometry and motion. This study images the fault strands that likely ruptured in the mainshock with higher spatial resolution than previous studies of the aftershocks.

The 2011 Mineral, Virginia mainshock was a complex rupture comprised of multiple subevents in a small fault area centered at depths between approximately 6.0 km and 8.0 km [2, 28, 29]. Previous studies found that the aftershocks occurred mostly within a tabular zone exhibiting strike and dip similar to one of the nodal planes of the mainshock focal mechanism [2, 16, 20, 24, 30]. Previous studies also identified aftershock locations approximately 10–20 km to the east and northeast of the mainshock that were interpreted to be triggered by stress transfer due to the mainshock [20, 24, 31]. The previously determined regional moment-tensor solutions for 16 of the largest aftershocks show a diversity of nodal plane orientations, although most of the focal mechanisms depict thrust motion [16, 20].

Shortly after the August 23, 2011, mainshock, teams from Virginia Tech, University of Memphis, Lehigh University, Incorporated Research Institutions for Seismology (IRIS) and the United States Geological Survey (USGS) rushed to the epicentral area and deployed a number of portable seismic instruments [20]. This aftershock survey recorded a rich dataset of aftershocks that is archived at the IRIS Data Management Center (http://www.iris.edu/dms/dmc/).

From August 23, 2011, through 31 December 2011 just over 1600 aftershocks were located by Dr. Martin Chapman, Qimin Wu, and Jacob Beale of Virginia Tech when the seismicity rate was high and the location capability was maximum and stable [30]. Chapman and Wu applied a detection algorithm that combined the short-term average and long-term average ratio (STA/LTA) method [32, 33] with a cross-correlation method that scanned through the continuous waveform recordings with signal templates to detect as many events as possible. The cross-correlation method refined the STA/LTA detections by removing many false triggers and detected additional small events. The automatic detections were visually reviewed and a total of 1600 aftershocks were detected within the 129-day time range. These 1600 well-recorded events were used as the dataset for this study.

Waveform data from weak-motion sensors collected by the Virginia Tech XY, IRIS Ramp YS, and USGS SY deployments of aftershock monitoring seismic stations were used in this study (Figure 1); data from strong-motion accelerometers deployed as part of the aftershock survey were not used in this study. A two-layer crustal velocity model developed for the CVSZ [2, 34] was used with both methods employed in this study.

Two separate analytical methods were utilized during this study to calculate with high spatial resolution the geometry of the fault structures illuminated by the aftershocks and determine the postmainshock fault deformation. The first analytical method, from Ebel et al. [35], is a modified version of the double-difference earthquake location scheme of Waldhauser and Ellsworth [36] and was used to calculate precise relative aftershock locations. The second analytical method is the moment-tensor inversion technique of Ebel and Bonjer [37], which was used to calculate focal mechanisms for dozens of individual aftershocks on the fault structures. Fault strikes and dips from both analytical methods were compared with each other to determine if the geometry of the fault structures depicted by both methods were consistent and to infer details regarding postmainshock fault deformation of the imaged fault structures. The calculated geometry of the fault structures was compared with known faults and preexisting planes of weakness near the location of the projected intersection of the fault structure with Earth’s surface to determine if the fault structures activated by the August 23, 2011, main event can be found.

A number of earthquakes from the aftershock sequence (called here “secondary” events) were located relative to spatially and temporally well-constrained “master” events using the technique outlined in Ebel et al. [35]. In the Ebel et al. [35] technique, cross-correlation was used to determine precise arrival time differences of P and S waves for the master event relative to a secondary event at common stations surrounding the epicentral area. These differences in the P and S travel times were used to calculate a highly precise offset in hypocentral location and origin time of the secondary event relative to the master event. Relative hypocentral differences less than 100 m are resolvable using this analytical method [35, 38]. To quantify the full range of the uncertainties in the relative hypocentral locations, a jackknife analysis was carried out. The jackknife method involves resampling a dataset numerous times while systematically omitting one individual observation per resample in order to determine the influence of each observation on the dataset as a whole. The full set of jackknife locations was used to calculate the standard deviation of the latitude, longitude, and focal depth of each relative hypocentral determination.

Several assumptions are made in the application of the relative location method of Ebel et al. [35]. The master and secondary events are assumed to have similar focal mechanisms, implying that the P and S waves that arrive at each station for the two events take off from similar locations on the body-wave radiation patterns from each source. Raypaths from the master event and secondary event to a common station are assumed to be similar, and thus, the waveforms from both events would experience the same scattering and frequency attenuation when recorded at a given station. The velocity structure in the hypocentral region that encompasses both the master and secondary events is assumed to be uniform, and any minor deviations in the ray paths to a given station for the two events are assumed to cause only small differences in the waveforms of the two events given the station sampling frequency of 100 Hz or 200 Hz for the data analyzed in this study.

Absolute aftershock locations of 1600 events provided by Virginia Tech were used to create multiple “clusters” that each contained one master event and tens to hundreds of possible secondary events for use in the relative location analysis (Figure 2), and a separate relative location analysis was carried out for the events in each cluster. The purpose of this was to maximize the chances that high-precision relative locations could be computed. An event near the center of each cluster that had a high signal-to-noise ratio and was well-recorded on as many stations as possible was chosen as the master event. The secondary events in a cluster were individually analyzed with the master event to determine if an accurate location of each secondary event relative to the master event could be calculated. For each master event-secondary event pair that was analyzed, a Butterworth filter (2-Hz and 10-Hz passband) was used to filter the signals before processing in the relative location analysis. Time windows encompassing the arrival times of the P waves (0.1 second before to 0.3 seconds after the hand-picked P-wave arrival) on the vertical component and S waves (0.5 seconds before to 1.5 seconds after the hand-picked S-wave arrival) on all three components were created for the master event and the secondary event. These time windows were used in the cross-correlation to find the precise difference in P- or S-wave arrival time that yielded the maximum similarity between the master event and secondary event P and S waveforms at a common station (Figure 3). Normalized cross-correlation coefficients (C-values) greater than 0.6 were assumed to represent cross-correlations with a high degree of similarity between waveforms, and therefore, cross-correlations of body waves at a common station with C-values less than 0.6 were not used to calculate relative locations (Figure 3). Of the three station components for which S-wave cross-correlations were performed (i.e. north-south, east-west, and vertical), the component with the largest C-value > 0.6 was used to calculate the S-wave arrival-time difference between the master event and the secondary event.

Some secondary event waveforms showed similar shapes but with opposite polarities when compared with the master event waveforms, resulting in large negative C-values determined in the cross-correlations. This can happen, for example, when two events are very close to each other but have opposite focal mechanisms (such as one thrust and one normal event, with both events having similar fault plane and slip orientations). For this reason, the maximum absolute C-values were used to determine which waveforms to include in the analysis (i.e. |C-value| > 0.6) as well as the lag time of the secondary event waveform relative to the master event waveform. Waveforms with |C-value| > −0.6, such as those shown in Figures 3(a) and 3(b) and Figure 3(c), were included in the relative location computations. In some cases, the waveforms for the master event, secondary event, or both were very sinusoidally shaped, which resulted in many cross-correlation peaks of similar values at different lag times. For these cases, it was unclear whether the maximum C-value represented the correct shift between the master event and secondary event even though the maximum C-value was greater than our acceptance value of 0.6. Under these circumstances, the time differences between the master events and the second events were not used for the relative location analysis due to this uncertainty in the correct lag time.

The lag times associated with the maximum C-values for body-wave cross-correlations at each station surrounding the events that passed the acceptance criteria described here were used to calculate the relative arrival time differences for the relative location analysis. In many cases, the RMS errors between the predicted and observed relative arrival times were below about 0.02 seconds for the relative locations, which indicates that the maximum possible precision of these relative locations was achieved.

Relative locations for secondary events in a cluster around individual master events were computed, and then the relative locations of pairs of events from different clusters were used to find accurate relative locations of the different clusters relative to each other. By doing this effectively, the master event of every cluster was relocated relative to one well-defined master event using the relative location method. This produced a highly precise reconstruction of the relative locations of the analyzed events across the entire aftershock zone.

Focal mechanisms were determined for a number of aftershocks at various locations in the aftershock zone using the moment-tensor inversion method of Ebel and Bonjer [37]. First-motion amplitudes and polarities were read from the direct P waves at stations surrounding an event. The S waves were rotated so that the amplitudes from the transverse components of the direct S waves could be read (SV waves were past the critical angle at the surface and so could not be used in the analysis). The instrument gains were used to convert the P and SH amplitudes to ground-motion amplitudes. The ground-motion amplitude data were input into a linear least-squares inversion to calculate the five independent components of a traceless moment tensor that best predicts the observed amplitudes and polarities. The resulting moment tensor was decomposed into the largest possible double-couple moment tensor and a residual compensated linear vector dipole (CLVD) moment tensor. The source moment, strike, dip, and rake were calculated from the double-couple. An additional nonlinear least-squares inversion code was used to find seismic moment, strike, dip, and rake values that optimize the fit between the observed amplitudes and those calculated from the double-couple solution.

This method was only appropriate for small magnitude events due to the assumption that the source time function is so short in duration that it can be represented by a Dirac delta function [37]. Some other assumptions and limitations are part of the application of this moment-tensor inversion method [37]. The seismic velocity structure between all sources and receivers was assumed to be the same for all events, and anelastic attenuation was assumed to be so small that it could be ignored in the analysis. All stations were located at far-field distances from the source, and only direct body waves were used in the moment-tensor inversions to avoid distortion from reflected and refracted waveforms. The seismic network data and velocity model used for the relative location analysis were also used in the moment-tensor analysis. Focal mechanisms were determined for events that had been relocated using the relative location method so that nodal plane orientations of the focal mechanisms could be compared with trends in the seismicity near the aftershock locations.

Several independent tests were conducted to determine the accuracy of the focal mechanisms found using the moment-tensor inversion method. Initially, inversions using only P-wave first-motion amplitudes were performed because the initial P-wave amplitudes were relatively easy to identify on many of the station seismograms. S-wave first arrivals can be masked by scattered P-waves and P-to-S wave conversions, making a clear S-wave first arrival difficult to determine [37]. Thus, as a first step in the focal mechanism procedure, P-wave-only inversions were performed using first-motion amplitudes read from common stations for events with small offsets in relative hypocentral locations, and the resulting focal mechanisms were compared to determine if the different events had similar double-couple solutions. S-wave first-motion amplitudes at common stations were then included in the inversion to determine whether the focal mechanisms found using both P- and S-wave amplitudes were similar to the focal-mechanism solutions using only P waves for the same events (Figure 4). S-wave first arrivals at stations closer than 3 km from the hypocenter were not used because the uncertainty in the hypocentral location can lead to large uncertainties in the rotation of the observed horizontal components to SV and SH ground motions. Also, small shifts in the hypocentral location could result in large changes in station take-off angle on the SH-wave radiation pattern.

In many cases, the addition of the S-wave amplitudes caused little change in the focal mechanisms relative to the P-wave-only inversion, or they caused a rotation of the nodal planes compared to the P-only solution that was consistent among events with nearby hypocenters (Figure 4). Additional body-wave amplitude picks at stations that were not common among the events with small offsets in hypocentral locations were then included in the inversion to further test the calculated focal mechanisms. The most accurate focal mechanisms were those that had little change in inversion results when independent data were included with common P- and S-wave stations, for which the RMS error between calculated and observed amplitudes was small, and for which the CLVD component was less than 1/6 of the double-couple component of the moment tensor.

Accurate focal mechanism solutions could not be determined for every event analyzed using the moment-tensor inversion method. An antialiasing filter used in the analog-to-digital conversion also appeared to create acausal noise before the P-wave [39, 40], making the first motion difficult to determine for some stations, especially for low-amplitude P-wave first arrivals (Figure 5). Also, the radiation patterns of events east of the main fault could not be constrained due to the relatively limited azimuthal coverage of the focal spheres due to the locations of the events relative to the recording stations.

The resulting relative location hypocenters reveal trends in seismicity that were not apparent from the absolute hypocenter locations alone, illuminating several major seismicity lineations within the main zone of aftershocks (Figure 6). These trends are interpreted in Figure 7 as two larger fault planes (Q1 and Q2) and one smaller fault plane (Q3), and for each of these, a best-fit plane was computed through the trend of the seismicity that is interpreted to lie along an individual fault plane (Figures 7 and 8). Surface projections of the two larger fault planes (Q1 and Q2) are shown on Figure 9 in addition to the surface projection of a single best-fit plane that was computed by using all of the aftershocks in the primary aftershock cluster (Q). Off-fault hypocenters near the intersection of the primary fault planes were not used in the calculation of the geometry of the primary fault planes in (Figure 8) to ensure that accurate strikes, dips, and surface projections of each interpreted plane could be determined.

Fifty-nine focal mechanisms were determined using the moment-tensor inversion method for events that had been relocated using the relative location method (Figure 10). Accurate focal mechanisms were created for events with hypocentral depths ranging from 2.45 km to 8.00 km and moment magnitudes (M_W) ranging from 0.48 to 2.69. Uncertainties in the calculated strikes and rakes of the focal mechanisms in this study range from ±0.01° to ±6.35°, whereas uncertainties in the calculated dips range from ±0.01° to ±1.61°. Focal mechanisms calculated for events on the two primary fault planes of the main fault structure generally depict thrust fault motion and have one nodal plane that corresponds closely with a fault plane delineated by the relative location hypocenters (Figure 10). Instances of normal-fault motion on planes with similar strikes and dips as nearby thrust focal mechanisms were calculated for two events on the main fault structure (Figure 10). In addition, normal faulting also occurred for several events that were located about 1–2 km northwest of the main aftershock zone (Figure 10). These events are located in the footwall block northwest of the main fault, and they seem to occur on a shallow planar feature and have similar strikes and dips as nearby thrust events on the main fault.

For the 2011 Mineral, Virginia aftershock zone, previous studies have interpreted the fault structure as one or two planar features (e.g. [24, 41]); however, in this study, it is determined that the aftershock fault structure is comprised of the intersection of two major fault structures along with probably a third smaller structure. The relative location analysis finds a complicated set of intersecting fault structures that consist of two primary planes (Q1 and Q2), and one possible smaller planar feature (Q3), along with some off-fault seismicity near the intersection of planes Q1 and Q2 (Figure 8). Fault planes Q1 and Q2 encompass most of the primary aftershock zone, albeit with noticeably different fault strikes. Smaller fault plane Q3 is interpreted because it best explains the relative location pattern in the southwestern part of the primary aftershock zone, a pattern that has a fault strike close to that of plane Q1. Planes Q1, Q2, and Q3 are all part of the Quail Fault Zone named by Horton et al. [17, 18, 24].

The best-fit planes calculated for the each of the three fault orientations imaged by the relative location results produced a strike and dip for each of the fault planes. The best-fit plane for the northernmost primary fault plane (Q1) has a strike and dip of 045°/67° SE, whereas the best-fit plane for the southernmost primary fault plane (Q2) has a strike and dip of 002°/72° SE. The strike and dip of Q3 is 056°/58° SE as determined by the trend of the relative aftershock locations, striking northeast similar to Q1 but with a shallower dipping plane.

Results of the moment-tensor inversion analysis for events on the northern primary fault plane (Q1) and southern primary fault plane (Q2) show focal mechanisms with one of their nodal planes with similar strikes and dips to the trends in seismicity as imaged in the relative location analysis (Figure 10). The focal mechanisms for events that were interpreted to be associated with the two primary fault planes were averaged to evaluate if the geometry of the fault structures depicted by both methods were consistent. The average strike, dip, and rake of the probable fault plane of the focal mechanisms calculated from the mean poles to the fault planes for events associated with Q1 are 42.99° ± 2.04° / 57.10° ± 0.73° SE / 90.15° ± 1.56° and the average strike, dip, and rake of the probable fault plane of the focal mechanisms calculated for events associated with Q2 are 30.88° ± 1.48° / 73.02° ± 0.24° SE / 93.23° ± 1.39°. No focal mechanisms were calculated for events associated with plane Q3.

The strike and dips calculated for the northernmost primary fault plane (Q1) using the relative location and moment-tensor inversion results generally depict similar trends: the strike of Q1 was calculated as 043° and 045° and the dip of Q1 was calculated as 57° SE and 67° SE. The dips calculated for the southernmost primary fault plane (Q2) are also similar and were calculated as 72° SE and 73° SE; however, the strikes calculated for Q2 were not as similar. The strike of the best-fit plane calculated using the relative aftershock locations (002°) was more northerly than the average strike of the focal mechanisms (031°) that were calculated for events on the same fault plane. The northeast striking and southeast dipping nodal planes of the focal mechanisms for events on primary fault plane Q2 consistently trend northeast rather than north as expected by the strike of the best-fit plane calculated for Q2.

The two different strikes calculated for Q2 are interpreted to depict a complex fault structure that resulted in structural control of the aftershock distribution for the events associated with Q2. The northern-striking linear trend in seismicity depicted by the relative aftershock locations (Q2) is interpreted to represent a highly fractured geologic unit that is located between two relatively unfractured geologic units, the westernmost potentially being the Ellisville Pluton. The contacts between the geologic units are interpreted to strike northward and seismicity was confined to the highly fractured north-south striking geologic unit. We interpret the joints or fractures in the highly fractured geologic unit to be oriented approximately parallel to one another and strike northeast as inferred from the strikes of the probable fault planes of the numerous focal mechanisms calculated for events associated with Q2. These parallel northeast-striking joints experienced slip during the aftershock sequence following the August 23, 2011, mainshock. Seismicity has been shown to be confined to fractured geologic formations in other seismically active regions [42], and focal mechanisms for aftershocks have been calculated with consistent but different nodal plane orientations than the nodal planes calculated for the mainshock for aftershock sequences in other seismic zones [43]. The combination of the relative location and moment-tensor inversion results allows for high-resolution interpretations of complex fault kinematics and structure.

McNamara et al. [20] interpreted the aftershock hypocenters as occurring on a single-fault plane with a strike and dip of 036° ± 12° / 050° ± 6° SE and Horton et al. [24] calculated a single-fault plane with a strike and dip of approximately 035°/48° SE. Horton et al. [24] also noted that the fault geometry of the aftershock zone appeared to be “concave-east” as a result of the intersection of two planar features with different orientations and calculated strikes and dips for each of the two planar features. Horton et al. [24] calculated a strike and dip of approximately 047°/59° SE for the northern part of the fault structure (similar to Q1 in this study) and a strike and dip of approximately 029°/62° SE for the southern part of the fault structure (similar to Q2 in this study). Horton et al. [24] did not identify a third planar feature (Q3 in this study); therefore, the southwest planar feature in Horton et al. [24] likely represents a combination of Q2 and Q3. The strike and dips of the two planar features identified by Horton et al. [24] are similar to the orientations of Q1 and Q2, identified in this study, with the exception of the strike of the southernmost planar feature (Q2) that was calculated using the best-fit plane of the relative location results. The average strike of the probable fault planes of the focal mechanisms associated with Q2 was 031°, similar to the strike of 029° calculated by Horton et al. [24].

The shallow cluster of events offset about 1–2 km to the northwest of the main fault illuminate a small planar fracture with strike and dip that are similar to those for the main fault, but this northwestern seismicity cluster seems to have experienced normal fault motion in the month following the mainshock (Figure 10(c)). Coulomb stress transfer analyses of the mainshock calculated positive Coulomb stress change at the location and depth of the northwestern seismicity cluster, indicating stress increase and reverse faulting would be expected at this location [30, 44]. Differences between the orientations of the receiver faults used in the Coulomb stress transfer analyses and the orientation of the nodal planes of the focal mechanisms in the northwestern seismicity cluster that were calculated in this study may explain how normal faulting could occur within the calculated location of positive Coulomb stress change. The strike and dip of the receiver faults used in these calculations were 026°/55° SE [44] and 340°/60° NE for receiver faults located within 0 to 4 km below ground surface [30], the same depth interval of the northwestern seismicity cluster. In this study, the orientation of the southeast dipping nodal planes from the four focal mechanisms calculated for aftershocks on the northwestern seismicity cluster range in strike from approximately 033° to 035° and range in dip from approximately 18° to 47° SE with an average strike and dip of 034°/33° SE. This normal fault motion, somewhat offset from the main fault structure, could be due to brittle rupture of a region in the mainshock footwall that experienced tensional stresses following the main thrust event [45] (Figure 11).

Relative location hypocenters of the aftershocks several kilometers to the northeast of the main fault depict a steeply dipping planar feature with a strike and dip of 033°/85° SE, and this structure has been named the Fredericks Hall Fault by Horton et al. [17, 18] (Figure 9). The radiation patterns of the body waves could not be constrained by the moment-tensor inversion method for events associated with the Fredericks Hall Fault due to the lack of stations east and north of the hypocenters. The hypocentral locations of the small cluster of events between the main fault aftershock zone and the Fredericks Hall Fault show little in the way of a consistent trend in strike or dip of a planar feature. For mainshocks at other continental locations outside of the CVSZ, distributed aftershock zones and abundant off-fault seismicity similar to the events northeast of the main fault and the shallow offset cluster of events to the northwest following the 2011 main event in Virginia have been observed for blind thrust rupture events [45] and are likely a result of redistribution of stress following the mainshock. Thus, the off-fault clusters of aftershocks in Virginia following the 2011 earthquake are probably not unusual for the aftershocks of intraplate earthquakes in continental settings.

The location of the Mineral, Virginia, earthquake within the CVSZ may have resulted from a combination of both the east-west to southeast-northwest direction of maximum compressive stress [2, 15, 16, 20] and the northeast striking geophysical anomalies near the epicentral area [24, 27]. The increase in the shear stress-to-normal stress ratio due to the bend in strike of the geologic units near the main event and the southeast-northwest orientation of maximum compressive stress near the mainshock epicenter likely controlled the location of fault slip [27]. Focal mechanisms calculated for aftershocks in this study depict thrust motion with P axes ranging from east-west to southeast-northwest, indicating a direction of average maximum compressive stress similar to that found in previous studies [2, 15, 16, 20] (Figure 10). The primary source of the east-west and southeast-northwest regional compressive stress in the CVSZ has been attributed to the opposition of ridge-push forces due to cooling of the lithosphere created at the Mid-Atlantic Ridge and the continental drag forces created by the basal drag friction between the continental lithosphere and the asthenosphere below the North American Plate [46].

The geologic units of the CVSZ have experienced multiple metamorphic and tectonic events during their geological history, which resulted in foliations, structural fabric, faults, and dikes that all contribute to preexisting orientations of structures in the epicentral region [19, 22, 23, 25-27]. North-to-northeast striking foliations and joints have been found at the surface [19, 22, 25] and preexisting faults have been identified in the subsurface [9, 26]. Dips of faults observed on the Interstate-64 (I-64) seismic profile are shallower than the dips determined for the Quail Fault structure, indicating that the Quail Fault does not extend as far south as I-64 as a structure with visible offset on the seismic profile [26]. Rock outcrops between the surface projection of the Quail Fault and the southern extension of the Long Branch Fault have structural fabrics that are indicative of a past high-strain environment, and this area has been labeled the Bend of River high-strain zone [25] (Figure 9). Hughes et al. [25] identified two distinct planar foliation orientations in the Bend of River high-strain zone with strikes and dips of 001°/25° SE and 037°/42° SE. The latter has been observed to be parallel to the east and southeast-dipping, compositionally variable layers and discontinuous lenses comprising the Chopawamsic Formation [25]. The strikes of these foliations are similar to the north- and northeastern-oriented strikes of 043° to 045° for primary fault plane Q1 and 002° to 031° for primary fault plane Q2 as determined in this study. Joints in the Chopawamsic Formation overlying the epicentral area have been found to strongly trend N15°E and N50°E to N60°E and to dip southeast [22]. Although no evidence of shear slip along these joints has been identified, some evidence of small-scale slip has been observed on numerous foliation-joint planes with similar strikes and dips as the two primary planes of the Quail Fault structure determined in this study [22].

North to northeast trending magnetic lineaments identified by an airborne geophysical survey flown over the epicentral region of the Mineral, Virginia, earthquake were interpreted to be folded and faulted Paleozoic metavolcanic and metasedimentary formations [27]. Shah et al. [27] also interpreted numerous north to northwest trending magnetic anomalies that cross-cut the Paleozoic era lineaments as Jurassic dikes. The Jurassic dikes and Paleozoic formations share similar strikes as primary fault planes Q2 and Q1, respectively. The contact between the dikes and adjacent geologic units are considered a plane of weakness, and dike-parallel faults and joints are often observed near dike contacts [47]. In support of these observations, mining-induced earthquakes referred to as “rock bursts” are known to occur at the interface between dikes and the adjacent geologic units [48]. It is likely that the geometry of the faulting on the Quail Fault and the fault structures activated during the aftershock sequence was controlled by a preexisting structural fabric that was activated during the 2011 Mineral, Virginia, mainshock.

The results of this study add to the growing body of evidence that earthquake ruptures in central and eastern North America (CENA) can be complex in space, time, or both. The M_w 5.1 August 2020 earthquake at Sparta, North Carolina showed a complex fault structure with multiple fault planes and focal mechanisms along with surface faulting [49]. Older earthquakes in CENA, such as the 1982 M5.8 Miramichi, New Brunswick earthquake and the 1988 M5.9 Saguenay, Quebec earthquake along other events in the region, were found to have heterogenous ruptures in time and space (e.g. [50, 51]). The results found in this study of the aftershocks of the 2011 Mineral, Virginia, earthquake show that high-resolution studies of the relative locations and focal mechanisms of an aftershock sequence can reveal details of the mainshock rupture that are not found by more standard analyses (such as [24, 41]). If the strong ground motions from CENA earthquake are to be properly modeled, details of the earthquake rupture histories must be known. Furthermore, details of which fault segments ruptured in the mainshock, such as those inferred in this study, can provide a better understanding of what fault structures may become seismically active in CENA in the future.

The relative location method used in this study defines a complicated fault structure that was illuminated by aftershocks that followed the Mineral, Virginia, 2011 main event. The moment-tensor analysis in this study produced focal mechanisms that have nodal planes with similar strikes and dips as the trends delineated by the seismicity. The majority of the focal mechanisms depicts thrust faulting with a few instances of normal faulting on the main fault structure and normal faulting on the shallow offset structure slightly northwest of the main fault structure. Results of the relative location method and the moment-tensor inversion method suggest that the active fault structure is comprised of three fault planes, two primary planes and one smaller plane with strikes and dips ranging from 043° to 045° / 57° to 67° SE for Q1, 002° to 031° / 72° to 73° SE for Q2, and 056°/58° SE for Q3.

Trends in the structural fabric of the bedrock in the epicentral area have similar strikes and dips as those calculated for the two primary fault planes and for the smaller fault plane on the main fault structure. The preexisting structural fabric within the bedrock of the CVSZ, and possibly the intersection of north-trending Jurassic dikes and northeast-trending Paleozoic geological contacts, likely controlled the geometry of the multiplanar Quail Fault Zone that was activated during the 2011 Mineral, VA mainshock. These preexisting weaknesses in the crystalline bedrock of the CVSZ likely control the locations of stress concentrations and the location of the initiation of large fault ruptures such as that in the August 23, 2011, Mineral, Virginia, earthquake.

The results of this study give an example of the type of complicated fault structures that can be activated during seismic events in the CVSZ and perhaps could be activated at other intraplate settings. The rapid deployment of dense seismometer networks and in-depth analyses of the recorded aftershock sequences following mainshocks as exemplified in this study are crucial for the advancement of the scientific understanding of active structures associated with the infrequent, yet damaging intraplate earthquakes in the central and eastern U.S.

The absolute locations of the aftershocks, P and S arrival time picks and aftershock waveforms from the portable seismic instrumentation were obtained from the seismology group at Virginia Tech. Additional waveforms for the relative location analysis were obtained from the SAGE Data Management Center operated by the EarthScope Consortium. Aftershock relative locations and focal mechanisms that were determined in this study are available from the authors

The authors declare no conflicts of interest regarding this work.

The material in this paper is based upon work supported by the U.S. Geological Survey under Grants Nos. G13AP00043 and G13AP00044.

We thank Dr. Martin Chapman and Qimin Wu of Virginia Tech for providing us with the waveform data, P- and S-wave arrival time picks, and the preliminary aftershock locations. We thank Andrea Billi and two anonymous reviewers for their comments on an earlier draft of this article. The contents of this article are based on the M.S. thesis research of Stephen Hilfiker in the Department of Earth and Environmental Sciences at Boston College. That thesis, which is entitled High-Resolution Spatial and Temporal Analysis of the Aftershock Sequence of the August 23, 2011, M_w 5.8 Mineral, Virginia Earthquake, can be found at the persistent link http://hdl.handle.net/2345/bc-ir:107179 at Boston College.