Studying the impact of strong near-field and far-field earthquakes on earthfill reservoir dams: an integrated satellite and instrumental monitoring approach

2.1. Seismic intensity attenuation with epicentral distance

The observed intensity $I$ at a site located at hypocentral distance R from an earthquake of epicentral intensity $I_0$ can be estimated with the macroseismic attenuation relationship of the classical Kövesligethy-Blake type [17]:

where: $I\left(R\right)$ – seismic intensity on the Medvedev-Sponheuer-Kárník macroseismic scale (MSK-64, 1964), a 12-degree ordinal intensity scale [16]; the value denotes the intensity degree (“point”) and is not a 64-point count; $I_0$ – epicentral intensity [MSK-64 degrees]; for the November 3, 2025 Afghanistan earthquake [7, 15] $I_0=$ 7-8; $R$ – hypocentral (slant) distance from the focus to the site [km], computed as $R=\sqrt{{\mathrm{\Delta }}^2+h^2}$; $\mathrm{\Delta }$ – epicentral distance [km]; for the Surkhan dam site $\mathrm{\Delta }\approx$ 146 km; $h$ – focal (hypocentral) depth [km]; $h=$ 30 km for the November 3, 2025 event; $\alpha$ – anelastic (inelastic) attenuation coefficient [km^-1]; for the Surkhan region $\alpha =$ 0.004 km^-1 (a value adopted in this study in the absence of a published regional calibration; its influence on the predicted intensity is examined in the sensitivity analysis, Table 5).

The first term $-3·{\mathrm{l}\mathrm{o}\mathrm{g}}_{10}(R/h)$ accounts for geometric (spherical) spreading of seismic energy, which decreases as the wave front expands. The second term $-3\alpha (R-h)/\left(2\mathrm{l}\mathrm{n}10\right)$ describes the additional energy loss due to anelastic absorption in the crustal material. Substituting the event parameters yields $I\approx$ 5 MSK-64 at the Surkhan dam site, consistent with instrumental records (Fig. 2).

Fig. 1General flowchart of the analysis methodology for seismic-safety assessment of the reservoir dam

2.2. Ground motion prediction equation (GMPE)

The Peak Ground Acceleration (PGA) at the dam site is estimated using a generalized illustrative attenuation relation of the GMPE type. Eq. (2) shows, for background, the generalized functional form of empirical ground-motion models of the Boore-Atkinson family [9]. The quantitative prediction used in this study is, however, computed with the fully specified NGA-West2 model of Boore, Stewart, Seyhan and Atkinson (BSSA14) [19], implemented in the open-source pygmm library (version 0.8.0, global region setting, unspecified fault mechanism, default basin-depth handling; the model output is the median of the orientation-averaged horizontal component): for $M_w=$ 6.2 [15], Joyner-Boore distance $R_{jb}\approx$ 146 km (epicentral distance used as the far-field proxy) and VS30 = 350 m/s the median PGA is 0.012 g (Fig. 3):

where: PGA – peak ground acceleration [g]; recorded values at the Surkhan dam: horizontal 0.17 g, vertical 0.18 g; $M$ – earthquake magnitude; $M=$ 6.8 is the value reported by the national seismic network for the November 3, 2025 Afghanistan event (magnitude type not specified in the platform record; the USGS moment magnitude is $M_w$ 6.2 [15], and the automatic SeisComP solution of the national network displays MLv 7.0, Fig. 14); $R_{rup}$ – closest distance from the site to the rupture plane [km]; since the rupture geometry is not available, the hypocentral distance $R$$=\sqrt{{\mathrm{\Delta }}^2+h^2}\approx \sqrt{{146}^2+{30}^2}$ ≈ 149 km is used here as a far-field proxy for $R_{rup}$; $V_{S30}$ — time-averaged shear-wave velocity in the top 30 m of soil at the site [m/s]; $V_{ref}$ – reference rock velocity; $V_{ref}=$ 760 m/s (NEHRP class B/C boundary); $c_1$-$c_2$ – empirical regression coefficients of the generalised form; Eq. (2) is shown as schematic background only, and no numerical result in this paper is derived from it – the quantitative prediction is computed with the fully specified BSSA14 model [19] as stated above; $\epsilon$ – residual error term representing record-to-record variability; $\epsilon ~N\left(0,{\sigma }^2\right)$.

Fig. 2Seismic intensity IΔ as a function of epicentral distance for the November 3, 2025 Afghanistan earthquake (h= 30 km, α= 0.004 km-1). The Surkhan dam site (Δ≈ 146 km) falls at I≈ 5 MSK-64 for I0 = 7

Fig. 3Median PGA versus Joyner–Boore distance from the fully specified BSSA14 ground-motion model [19] (VS30 = 350 m/s, unspecified mechanism; computed with the pygmm v0.8.0 implementation) for Mw 5.5-7.5. The event point (Mw= 6.2, Rjb≈ 146 km, median PGA = 0.012 g) is marked; the recorded PGA (0.17 g) and the design limit (0.35 g) are shown as horizontal lines

The term $c_{2}M$ controls the source scaling – larger magnitude earthquakes radiate more energy across a broader fault area. The term $c₃\mathrm{l}\mathrm{n}(R_{rup} +...)$ captures the combined effect of geometric spreading and path-dependent attenuation. The last coefficient $c₆\mathrm{l}\mathrm{n}(V_{S30} /V_{ref})$ introduces a site-amplification correction: softer soils (lower $V_{S30}$) amplify ground motion relative to reference rock. For the Surkhan dam body (compacted clay-gravel fill, $V_{S30}$ ≈ 320-420 m/s), a modest amplification factor of approximately 1.2-1.4 is expected.

2.3. Equation of motion for the dam under seismic excitation

The dynamic response of the dam structure is governed by the classical multi-degree-of-freedom (MDOF) equation of motion (Fig. 4). When the dam is treated as a discretised finite-element system subjected to base acceleration ${\ddot{u}}_g\left(t\right)$ from the earthquake record, Newton’s second law yields [10] (Eq. (3)):

where: $\left[M\right]$ – global mass matrix [kg]; assembled from element mass matrices $m_e=\rho \int N^{T}NdV$, with $\rho$ the material density; $\left[C\right]$ – global damping matrix [N·s/m]; defined via Rayleigh proportional damping (see Eq. 4); $\left[K\right]$ – global stiffness matrix [N/m]; assembled from element stiffness matrices $k_e=\int B^{T}DBdV$, where $B$ is the strain-displacement matrix and $D$ is the elasticity tensor; $\left\{u\left(t\right)\right\}$ – nodal displacement vector [m]; $\left\{\dot{u}\left(t\right)\right\}$ = first time derivative (velocity), $\left\{\ddot{u}\left(t\right)\right\}$ = second derivative (acceleration); $\left\{r\right\}$ – influence (participation) vector; maps the scalar base acceleration to all structural degrees of freedom; $\ddot{u}g\left(t\right)$ – recorded ground acceleration time history [m/s²]; obtained from the 18 accelerogram channels at the dam.

Fig. 4Lumped-mass multi-degree-of-freedom (MDOF) idealisation of the dam, showing storey masses mi, inter-mass stiffness ki and damping ci, base excitation u¨gt, and nodal responses uit governed by Eq. (3)

Eq. (3) is solved numerically by direct time integration using the Newmark-$\beta$ method with parameters $\beta =$ 0.25 and $\gamma =$ 0.5 (constant average acceleration, unconditionally stable). The time step is chosen as $\mathrm{\Delta }t\le T_{min}$/10, where $T_{min}$ is the smallest period of interest in the structural response.

2.4. Finite element model – structural description

The boundary conditions of the model are a fixed base ($u=v=$0), lateral rollers ($u=$ 0) and a hydrostatic load applied to the upstream face at reservoir level 78 m. The discrete system of Eq. (3) is constructed by assembling a two-dimensional plane-strain finite element model of the Surkhan dam cross-section in the SAP2000 structural-analysis software. The main model parameters are summarised below. Mesh convergence was verified by progressively refining the discretisation until the change in the computed peak nodal displacement between successive refinements became negligible, ensuring that the reported response is independent of the mesh density (Fig. 5).

Fig. 5Surkhan dam FEM cross-section mesh (plane-strain, 240 elements, ~270 nodes, 540 DOF). Material zones, boundary conditions, hydrostatic loading (blue arrows), and seismic base excitation u¨g(t) are shown

Element and node count. The dam cross-section is discretised into 200 four-node isoparametric quadrilateral elements (quad4) and 40 three-node triangular elements (tri3) in transitional zones, giving a total of 240 elements and approximately 270 nodes. Each node carries two translational degrees of freedom (horizontal u and vertical v), resulting in $N=$ 2×270 = 540 active DOF in the assembled system matrices $\left[M\right]$, $\left[C\right]$, and $\left[K\right]$.

Material zones. The dam body is divided into four material zones: Zone 1 (clay core, $E=$ 40-60 MPa, $\nu =$ 0.40, $\rho =$ 1900 kg/m³); Zone 2 (filter/transition, $E=$ 60-90 MPa, $\nu =$ 0.35, $\rho =$ 1850 kg/m³); Zone 3 (gravel shell, $E=$ 80-120 MPa, $\nu =$ 0.30, $\rho =$ 1800 kg/m³); Zone 4 (rock foundation, $E=$ 500 MPa, $\nu =$ 0.25, $\rho =$ 2200 kg/m³). Linear elastic behaviour is assumed, consistent with the small-strain response at PGA = 0.17 g.

Boundary conditions. The base of the foundation block is fully fixed ($u=v=$ 0). Lateral boundaries use horizontal rollers ($u=$ 0, $v$ free). The upstream face carries triangularly distributed hydrostatic pressure for reservoir level HW = 78 m.

Seismic forcing scheme. Earthquake excitation is applied as uniform base acceleration ${\ddot{u}}_g\left(t\right)$ transmitted through the influence vector $\left\{r\right\}$ to all base-level DOF. Six accelerometers at the dam provided 18 accelerogram records on November 3, 2025; the dominant horizontal component (PGA = 0.17 g) is used. Time integration uses the Newmark-$\beta$ constant-average-acceleration scheme ($\beta =$ 0.25, $\gamma =$ 0.5, $\mathrm{\Delta }t=$ 0.01 s). The finite-element model is two-dimensional (plane strain), so the excitation is applied as a single horizontal component at the base rather than as a three-component motion. The trace plotted in the top panel of the following figure is a representative amplitude-modulated signal whose peak amplitude is scaled to the horizontal PGA of 0.17 g registered at the dam on November 3, 2025; it is used to demonstrate the stability and convergence of the time-integration scheme and does not reproduce the recorded waveform, which is shown separately in Fig. 37. Because this representative excitation sustains harmonic input over the whole record, the crest displacement it produces is a demonstration value of the order of tens of millimetres and is therefore much larger than the displacement of ±0.16 mm actually registered during the November 3, 2025 event (Fig. 6).

2.5. Rayleigh proportional damping

For earthfill dam materials, direct measurement of the damping matrix $\left[C\right]$ is not feasible. The Rayleigh proportional damping model expresses $\left[C\right]$ as a linear combination of $\left[M\right]$ and $\left[K\right]$:

where: ${\alpha }_R$ – mass-proportional Rayleigh coefficient [s⁻¹], controlling damping at low frequencies; ${\beta }_R$ – stiffness-proportional Rayleigh coefficient [s], controlling damping at high frequencies.

Fig. 6Newmark-β time integration of the dam dynamic equation of motion (β= 0.25, γ= 0.5, Δt= 0.01 s): two-dimensional plane-strain model. Top panel: input ground acceleration; middle panel: computed crest displacement ut(middle), and velocity u˙t (bottom). The input trace is a representative amplitude-modulated excitation scaled to the horizontal PGA of 0.17 g recorded at the dam on November 3, 2025; the instrumental record of the event itself is presented in Fig. 37. The dashed reference line marking the recorded maximum displacement (±0.16 mm) falls on the horizontal axis of the middle panel because that value is more than two orders of magnitude smaller than the vertical scale of the simulated response, which is governed by the sustained representative excitation

a) Input ground acceleration – representative synthetic trace scaled to the recorded PGA 0.17 g

b) Dam crest displacement response (demonstration of the integration scheme)

c) Dam crest velocity response

These coefficients are determined from the specified modal damping ratio $\xi$ at two natural frequencies ${\omega }_i$ and ${\omega }_j$:

where: $\xi$ – target modal damping ratio [dimensionless]; $\xi =$ 0.05 (5 %) for compacted earthfill dam materials under earthquake loading (ICOLD, 2016) [12]; ${\omega }_i$, ${\omega }_j$ – angular natural frequencies [rad/s] of the first and second dominant vibration modes of the dam (Fig. 7).

2.6. Stress-strain constitutive model for dam soil layers

The stress state within the dam body is related to strains through the linear elastic (Hooke’s) constitutive law in full three-dimensional tensor form:

where: ${\sigma }_{ij}$ – Cauchy stress tensor [Pa]; diagonal components ${\sigma }_{11}$, ${\sigma }_{22}$, ${\sigma }_{33}$ are normal stresses; off-diagonal ${\sigma }_{ij}$ ($i\ne j$) are shear stresses; ${\epsilon }_{kl}$ – infinitesimal strain tensor [dimensionless]; ${\epsilon }_{kl}=\frac{1}{2}\left(\frac{\partial u_k}{\partial x_l}+\frac{\partial u_l}{\partial x_k}\right)$; $C_{ijkl}$ – fourth-order elasticity (stiffness) tensor [Pa]; for an isotropic material it reduces to two independent parameters: Young’s modulus $E$ and Poisson’s ratio $\nu$.

Fig. 7Rayleigh proportional damping ratio ξω showing mass-proportional, stiffness-proportional, and total contributions. Target ξ= 5 % is exactly met at the two anchor frequencies ω1 and ω2 of the dam

For isotropic soil layers, Eq. (6) simplifies to the engineering matrix form $\sigma =D\epsilon$, where $D$ is the 6×6 elasticity matrix. In principal stress notation, the maximum shear stress is given by:

where: ${\sigma }_1$, ${\sigma }_3$ – major and minor principal stresses [Pa].

Fig. 8FEM stress distribution in the Surkhan dam cross-section (November 3, 2025, M= 6.8, PGA = 0.17 g) major principal stress σ1 (left) and maximum shear stress τmax=(σ₁-σ₃)/2 (right) [MPa]. The computed maxima are checked against the Mohr-Coulomb strength envelope in Fig. 9; that check is conditional on the assumed c and φ

The computed maximum stresses in the Surkhan dam body under the November 3 event were: major principal stress ${\sigma }_1=$ 1.63 MPa, minor principal stress ${\sigma }_3=$ 0.41 MPa, and maximum shear stress ${\tau }_{max}=(\sigma ₁-\sigma ₃)$/2 = 0.61 MPa, in agreement with the contour maxima of Fig. 8. These stresses are assessed against the strength of the fill in two ways. Relative to the nominal permissible stress adopted for the compacted fill in Table 1 they are an order of magnitude smaller. Relative to the Mohr-Coulomb strength envelope assumed for the fill ($c=$ 0.15 MPa, $\phi =$ 30°), however, the stress circle remains below the envelope by only about 0.03 MPa (Fig. 9); no plastic deformation is therefore predicted, but this margin is small and depends directly on the assumed strength parameters, which should be confirmed by site-specific laboratory and in-situ testing.

Fig. 9Mohr circle for the earthquake-induced stress state in the Surkhan dam body, plotted for the computed principal stresses σ1= 1.63 MPa and σ3= 0.41 MPa, which give τmax=(σ₁-σ₃)/2= 0.61 MPa in agreement with the contour maxima of Fig. 8. With the assumed Mohr–Coulomb strength parameters (c= 0.15 MPa, φ= 30°) the circle does not intersect the failure envelope, but the clearance is only about 0.03 MPa, so the mobilised shear is close to the assumed strength and the strength parameters require verification by site-specific testing

2.7. Seismic safety margin based on intensity

The seismic safety margin of the dam is expressed as the difference between the design-basis intensity and the intensity recorded during the event (Eq. (8)), $\mathrm{\Delta }I = I_{design}-I_{recorded}$. Because MSK-64 is an ordinal (degree-based) scale, a ratio of two intensity values has no quantitative meaning; the margin is therefore reported in intensity degrees rather than as a factor:

where: $I_{design}$ – design seismic intensity for the Surkhan dam [MSK-64 degrees]; $I_{design}$ = intensity degree 7, per the original design specifications; $I_{recorded}$ – recorded seismic intensity at the dam during the earthquake event [MSK-64 degrees]; $I_{recorded}$ ≈ intensity degree 5.

Substituting: $\mathrm{\Delta }I=$ 7-5 = 2 MSK-64 degrees. A margin $\mathrm{\Delta }I>$ 0 confirms that the intensity experienced at the site remained below the design basis; a margin of two full intensity degrees indicates that the loading was well below the design level. Quantitative capacity-to-demand ratios are meaningful only for the cardinal response quantities – acceleration, displacement and stress – and these are compared with their design limits in Fig. 10 and Table 1.

2.8. InSAR line-of-sight displacement decomposition

Satellite-based Interferometric Synthetic Aperture Radar (InSAR) measures surface displacement projected onto the satellite’s line-of-sight (LOS) direction. The LOS displacement d_LOS is related to the three-dimensional ground displacement components [8]:

where: $d_{LOS}$ – LOS displacement [mm]; positive toward the satellite (uplift or approach); $u_E$, $u_N$, $u_U$ – East, North, and vertical (Up) displacement components [mm]; $\theta$ – satellite incidence angle measured from vertical [degrees]; for Sentinel-1 over the Surkhan area $\theta \approx$ 33°-43°; $\phi$ – satellite heading (flight azimuth) angle [degrees]; approximately 349° for Sentinel-1 ascending orbit. A single line-of-sight geometry does not by itself allow the three displacement components $u_E$, $u_N$, $u_U$ to be separated; the horizontal and vertical series reported by the national platform result from the platform’s combination of acquisition geometries and processing assumptions, and are used here as provided.

Fig. 10a) Comparison of recorded response parameters with design limits and b) capacity-to-demand ratios for the cardinal response quantities – acceleration, displacement and shear stress. All ratios exceed 1.0, confirming adequate structural reserve

a)

b)

The interferometric phase $\mathrm{\Delta }\phi$ recorded by the satellite is directly proportional to $d_{LOS}$:

where: $\lambda$ – radar wavelength [m]; for Sentinel-1 C-band $\lambda =$ 0.0556 m (5.56 cm), giving sensitivity of $\lambda$/2 = 2.78 cm per 2$\pi$ phase cycle; $\mathrm{\Delta }\phi$ – interferometric phase difference [radians] between two SAR acquisitions.

The maximum LOS displacement magnitude observed at the Surkhan dam for the full monitoring period (January 2017-December 2023) was –79.6 mm (subsidence) and +86.8 mm (uplift), corresponding to phase changes of approximately –18.0 rad and +19.6 rad, respectively. Seasonal phase oscillations correlate directly with reservoir water level fluctuations, confirming the hydrostatic loading origin of the observed displacements (Fig. 11).

Table 1 summarises all measured seismic response parameters and the corresponding permissible limits specified in the design documentation of the Surkhan dam; the stress margin in the table is conditional on the assumed Mohr-Coulomb parameters $\left(c,\phi \right)$, which have not been measured on site.

Under the adopted parameter assumptions, the integrated application of Eqs. (1-10) shows that all measured response quantities remain within the permissible engineering limits of Table 1; the stress rows are conditional on the assumed Mohr–Coulomb strength parameters rather than on measured values. No plastic deformation, instability, or critical cracking was triggered by the event. The mathematical models presented here provide a quantitative framework for interpreting the combined satellite and instrumental monitoring data, and can serve as the basis for future probabilistic seismic hazard assessments of the dam.

Fig. 11a) InSAR line-of-sight geometry showing incidence angle θ≈ 38° and displacement component decomposition for Sentinel-1 ascending orbit; b) illustrative simulated LOS displacement time series reproducing the seasonal oscillation pattern of the instrumental InSAR record (which covers January 2017 – December 2023); the series is extended schematically to the November 2025 event for context and is not measured data beyond December 2023

a)

b)

2.9. Response spectrum analysis

The response spectrum provides the maximum seismic response (acceleration, velocity, or displacement) of a linear single-degree-of-freedom (SDOF) oscillator as a function of its natural period $T$, for a given damping ratio $\xi$. It is the fundamental tool for evaluating whether a structure's natural period falls in a critical excitation range. The spectral acceleration $S_a\left(T\right)$ is constructed using the Newmark-Hall amplification approach [13]:

where: $S_a\left(T\right)$ – spectral acceleration [g] for a period-T oscillator with damping $\xi$; $A_a$, $A_v$, $A_d$ – spectral amplification factors for acceleration, velocity, displacement regions; for $\xi =$ 5 %: $A_a=$ 3.21, $A_v=$ 2.31, $A_d=$ 1.82 (Newmark-Hall, 1982) [13]; $T_a$, $T_v$, $T_d$ – transition periods [s] separating the spectral regions; $T_a$ ≈ 0.1 s and $T_v$ ≈ 0.5 s for rock-like site conditions, and $T_d$ marks the long-period transition to the displacement-controlled region; for $T<T_a$ the spectrum rises from the PGA to the plateau

The spectral displacement $S_d$, which governs structural deformation demand, is obtained from Sa through the fundamental kinematic relationship:

where: $S_d\left(T\right)$ — spectral displacement [m]; with $S_a$ expressed in units of $g$, the conversion factor $g=$ 9.81 m/s² is applied, $S_d=g·S_a·{\left(\frac{T}{2\pi }\right)}^2$; equivalent to the maximum displacement of a SDOF oscillator

For the Surkhan earthfill dam the fundamental period adopted is $T_1\approx$ 0.45 s, taken from the modal analysis of the finite-element model ($f_1\approx$ 2.22 Hz, Section 4.3) [10]. The empirical cross-check of Eq. (14), $T_1=$ CT·H^3/4/√Es, evaluated with the adopted parameter ranges (CT ≈ 0.025-0.045, $H=$ 85 m, $E_s\approx$ 80-120 MPa), gives ≈ 0.06-0.14 s, well below the modal value; empirical formulas of this type are calibrated for other dam typologies, so they are used here only as an order-of-magnitude comparison, and the discrepancy is carried as a modelling uncertainty (Table 5). At $T_1=$ 0.45 s, the spectral acceleration reaches $S_a\left(T_1\right)$ ≈ 0.55 g and the spectral displacement $S_d\left(T_1\right)$ ≈ 27.5 mm. These values represent the maximum possible SDOF demand – the actual dam response is significantly smaller due to multi-mode behavior and 3D effects (Fig. 12).

Fig. 12a) Response spectrum SaT constructed for PGA = 0.17 g and ξ= 5 % damping, showing the dam fundamental period T1≈ 0.45 s in the spectral acceleration plateau; b) spectral displacement SdT=g·SaT·T/2π2, with g= 9.81 m/s2 converting Sa from units of g, indicating SdT1≈ 27.5 mm as the upper-bound displacement demand.

a)

b)

2.10. Natural frequency and modal analysis of the dam

The natural frequencies of the dam are determined from the undamped free-vibration eigenvalue problem derived from Eq. (3) with $\left[C\right]=$ 0 and no external loading:

where: ${\omega }_n$ – $n$-th circular natural frequency of the dam [rad/s]; $f_n={\omega }_n/2\pi$ [Hz], $T_n=$ 1/$f_n$ [s]; $\left\{{\phi }_n\right\}$ – $n$-th mode shape vector (eigenvector); describes the deformation pattern at frequency ${\omega }_n$.

For an earthfill dam of height $H$, the empirical formula for the fundamental period $T_1$ is:

where: $C_T$ – empirical coefficient that absorbs the unit conversion [s·MPa^1/2·m^-3/4]; values of 0.025-0.045 are reported for earthfill dams in the compilation of [12]; $H$ – dam height [m]; $H=$ 85 m for the Surkhan dam; $E_s$ – average shear modulus of dam material [MPa]; $E_s\approx$ 80-120 MPa for compacted clay-gravel.

The modal participation factor ${\mathrm{\Gamma }}_n$ measures how strongly each mode is excited by the base acceleration:

where: ${\mathrm{\Gamma }}_n$ – modal participation factor [dimensionless]; governs the amplitude of mode n in the total response. For earthfill dams, the first mode typically captures 75-85 % of the total seismic mass.

The maximum displacement $u_{max}$ at any point is computed by combining modal contributions using the Complete Quadratic Combination (CQC) rule:

where: ${\rho }_{ij}$ – cross-correlation coefficient between modes $i$ and $j$ [–]; depends on $\xi$ and the frequency ratio ${\omega }_j$/${\omega }_i$ (Der Kiureghian, 1981) [14]; $S_{d_i}$ – spectral displacement of mode $i$ from Eq. (12) [m].

For the Surkhan dam ($H=$ 85 m, $V_{S30}\approx$ 350 m/s), the estimated fundamental frequency is $f_1\approx$ 2.22 Hz ($T_1\approx$0.45 s) and $f_2\approx$6.67 Hz ($T_2\approx$0.15 s). These are consistent with published free-vibration records for comparable embankment dams [5] and place the fundamental period on the constant-acceleration plateau of the response spectrum ($T_a\le T_1\le T_v$).

2.11. Liquefaction potential assessment

Although the Surkhan dam body consists of compacted clay-gravel fill (not susceptible to liquefaction), the alluvial foundation materials in the river channel zones require assessment. The simplified procedure of Seed and Idriss (1971) [11] is applied, comparing the cyclic stress ratio (CSR) induced by the earthquake with the cyclic resistance ratio (CRR) of the foundation soil [11]:

where: CSR – cyclic stress ratio induced by the earthquake [dimensionless]; ${\sigma }_v$ – total vertical (overburden) stress at depth $z$ [kPa]; ${\sigma }_v={\gamma }_{soil}·z$, where ${\gamma }_{soil}$ is the unit weight of the soil [kN/m³], ${\sigma }_v^{\text{'}}$ – effective vertical stress at depth z [kPa]; ${\sigma }_v^{\text{'}}={\sigma }_v-u_w$, where $u_w$ is pore water pressure; $a_{max}$ – peak ground acceleration at the surface [m/s²]; $a\_\left(max\right) =$ 0.17 g = 1.67 m/s²; $r_d$ – stress reduction factor [–]; $r_d\approx$ 1.0-0.00765$z$ for $z\le$ 9.15 m, $r_d\approx$ 1.174-0.0267$z$ for 9.15 m $<z\le$ 23 m (Seed & Idriss, 1971) [11].

The cyclic resistance ratio CRR is estimated from standard penetration test (SPT) blow count $N_{160}$ corrected for overburden and energy efficiency:

where: CRR_7.5 – cyclic resistance ratio for a reference earthquake of $M=$ 7.5 [dimensionless]; $N_{160}$ – SPT blow count normalized to 100 kPa overburden and 60 % hammer energy efficiency.

For earthquakes of magnitude $M\ne$ 7.5, CRR is adjusted by the Magnitude Scaling Factor (MSF):

The factor of safety against liquefaction at each depth is then:

where: $F{S}_{liq}$ – liquefaction factor of safety [–]; $F{S}_{liq}>$ 1.0 indicates no liquefaction risk.

For the Surkhan dam foundation ($M=$ 6.8, $a_{max}=$ 0.17g, representative $N_{160}=$ 20-30 for medium-dense alluvial sand-gravel): MSF = 102.24 / 6.22.56 ≈ 1.63, evaluated with the USGS moment magnitude $M_w=$ 6.2 [15], since the magnitude-scaling factor of the NCEER procedure is defined in terms of moment magnitude [18]; CSR at $z=$ 5 m ≈ 0.65 × (90/72) × 0.17 × 0.96 ≈ 0.133; CRR_7.5 ≈ 0.29 (for $N_{160}=$ 25, evaluated from Eq. (18), the clean-sand base curve of the NCEER consensus procedure [18]); CRRM = CRR_7.5 × MSF ≈ 0.29 × 1.63 ≈ 0.47; ${FS}_{liq}$ ≈ 0.47 / 0.133 ≈ 3.6 ≫ 1.0 (the conservative alternative of applying the national $M=$ 6.8, whose magnitude type is not documented, would give MSF ≈ 1.28 and ${FS}_{liq}$ ≈ 2.8, likewise ≫ 1.0). This confirms that liquefaction of the foundation materials was not triggered by the November 3, 2025 earthquake event, consistent with the absence of post-earthquake settlements or lateral spreading observed at the site.

4.1. Model-observation comparison: predicted versus observed response

A central requirement for a credible monitoring-and-modeling framework is that the analytical predictions reproduce the instrumentally observed quantities. Table 3 compares the values predicted by the principal mathematical models (Eqs. (1, 2, and 12)) with the corresponding observed values recorded at the Surkhan dam during the November 3, 2025 earthquake.

The seismic-intensity attenuation model (Eq. (1)) predicts $I\approx$ 4.6 MSK-64 at the Surkhan site (epicentral distance $\mathrm{\Delta }\approx$ 146 km, $I_0=$ 7), in close agreement with the instrumentally recorded value of 5 MSK-64 (ratio 0.92). The ground-motion prediction is not obtained from the schematic Eq. (2) but from the fully specified BSSA14 model [19]: the median free-field PGA is 0.012 g for $M_w=$ 6.2 [15] at $R_{jb}\approx$ 146 km with VS30 = 350 m/s and unspecified mechanism – roughly fourteen times below the recorded 0.17 g (Fig. 3). Possible explanations for the gap between the model median and the recorded value include local site and crest amplification, the sensor position and processing chain, and event-specific source effects; these remain hypotheses until the raw station records are examined. The size of the gap is consistent with local site and structural (crest) amplification and/or event-specific source effects that the far-field free-field model does not represent; the recording conditions therefore warrant verification before a cause is assigned. The seismic-intensity model, by contrast, reproduces the observed MSK-64 intensity closely (ratio 0.92).

The numerical difference of about two orders of magnitude (a factor of roughly 170 – a numerical comparison only, not a validation ratio; Table 3 marks the two quantities as not directly comparable) between the spectral displacement $S_d\left(T_1\right)$ ≈ 27.5 mm (Eq. (12)) and the recorded crest displacement of 0.16 mm requires explicit clarification. These two quantities are not directly comparable. The spectral displacement Sd is the theoretical maximum response of a linear single-degree-of-freedom (SDOF) oscillator tuned to the dam fundamental period ($T_1\approx$0.45 s) and subjected to the full design-level response spectrum anchored at PGA = 0.17 g – it represents an upper-bound demand under sustained resonant excitation. The 0.16 mm value, by contrast, is the actual relative displacement recorded at the dam crest during the November 3 event. The reduction by a factor of about 170 from the SDOF upper bound to the observed response is attributable to several physical factors: (i) the far-field nature of the event ($\mathrm{\Delta }$ ≈ 146 km), which delivered only a small fraction of the design-level spectral energy at the dam period; (ii) strong radiation and material damping in the 85 m earthfill body ($\xi =$ 5 %), which suppresses resonant build-up; (iii) the three-dimensional, multi-modal response of the embankment, in which the first mode captures only 75-85 % of the seismic mass and higher modes partially cancel; and (iv) the fact that the recorded motion was a brief transient rather than a steady-state harmonic excitation. The spectral displacement should therefore be interpreted as a conservative design-envelope quantity, while the 0.16 mm value is the realized response; both are consistent with a dam responding well within its elastic range.

4.2. Relationship between long-term InSAR and seismic-transient displacements

A key interpretive point concerns the relationship between the two independent displacement datasets used in this study, which operate on fundamentally different time scales and capture different physical phenomena. The InSAR displacements (seasonal amplitudes up to ±86.8 mm vertical and ±110.7 mm horizontal over 2017-2023) are quasi-static deformations driven by the reservoir loading–unloading cycle: as the water level rises and falls seasonally (Fig. 35), the hydrostatic load on the dam body and foundation changes, producing predominantly seasonal deformation with superposed cumulative point-wise trends whose reversibility remains unverified with a characteristic period of months. The accelerometer displacements (0.16 mm horizontal, 0.15 mm vertical) are dynamic seismic transients produced by the passage of seismic waves during the ≈60 s of strong shaking on November 3, 2025, with a characteristic period of seconds.

These two regimes are complementary rather than contradictory. The InSAR record establishes the long-term pre-earthquake baseline (January 2017 – December 2023) and shows that the deformation is dominated by the seasonal reservoir loading-unloading cycle; the aggregate mean displacement is near zero (Fig. 19), while individual monitoring points exhibit cumulative trends of up to a few tens of millimetres over the seven-year record (Figs. 28-34), so point-wise trend assessment remains part of routine monitoring rather than a closed conclusion. The accelerometer record, superimposed on this baseline, demonstrates that even a strong regional earthquake adds only a sub-millimeter transient perturbation – roughly two to three orders of magnitude smaller than the routine seasonal movements the dam already accommodates safely. The integration of the two datasets thus provides a more complete safety picture than either could alone: the satellite data provide evidence of predominantly seasonal deformation behaviour during 2017-2023, while the instrumental data document the dynamic response to the recorded event. This combined interpretation is the principal methodological contribution of the present study.

4.3. Modal frequencies and comparison with the literature

The fundamental and second natural frequencies of the Surkhan dam, estimated from the modal analysis (Eq. (13-15)), are $f_1\approx$ 2.22 Hz ($T_1\approx$ 0.45 s) and $f_2\approx$ 6.67 Hz ($T_2\approx$ 0.15 s). Table 4 places these values in the context of measured and computed natural frequencies reported for comparable earthfill dams in the literature, confirming that the present estimates fall within the expected range for embankment dams of similar height and material composition.

4.4. Uncertainty and sensitivity analysis

The quantitative results depend on several input parameters that carry inherent uncertainty. A qualitative sensitivity analysis (Table 5) identifies how variations in the key parameters propagate to the predicted response quantities and to the intensity margin ΔI. Because the recorded response is far below the displacement and acceleration limits, the safety conclusion is insensitive to the plausible variation of the individual parameters considered; the exception is the material shear strength, where the Mohr-Coulomb clearance is only ≈ 0.03 MPa with unmeasured $c$ and $\phi$, so strength verification remains a priority.

Taking the most unfavorable combination of parameters (high attenuation reducing apparent intensity, soft-soil amplification raising PGA, low shear modulus lowering the fundamental frequency toward the spectral plateau, and reduced damping), the predicted demand quantities increase by at most a factor of about two, which still leaves the dam response one to three orders of magnitude below the design limits in Table 1. The intensity margin $\mathrm{\Delta }I=I_{design}-I_{recorded}$) = 7 − 5 = 2 degrees is governed directly by the recorded intensity and is therefore insensitive to modeling assumptions; even if the recorded intensity were one full MSK-64 degree higher ($I=$ 6), a positive margin of one degree would remain, preserving the safety conclusion.

4.5. Geotechnical model parameters and reference standards

Table 6 summarizes the geotechnical and material properties adopted for the finite-element model of the Surkhan dam, together with the boundary conditions and the reference standards used to assign site and material parameters.

Two international reference frameworks underpin the parameter selection. NEHRP (the U.S. National Earthquake Hazard Reduction Program site classification) defines site categories by the time-averaged shear-wave velocity in the upper 30 m ($V_{S30}$); the Surkhan foundation rock ($V_{S30}\approx$ 760 m/s) corresponds to the NEHRP class B/C boundary, which corresponds to the reference velocity $V_{ref}=$ 760 m/s appearing both in the generalised form of Eq. (2) and in the BSSA14 model [19]. ICOLD (the International Commission on Large Dams) provides the standard guidance for the dynamic analysis of embankment dams, including the recommended modal damping ratio ($\xi =$ 5 %) for compacted earthfill under earthquake loading and the empirical site coefficient $C_T$ for the fundamental-period formula (Eq. (14)). Adopting these widely accepted standards rather than locally calibrated values makes the methodology transferable to other dams and the results comparable with international practice.

Finally, the question of reservoir-triggered seismicity (RTS) deserves comment. The Surkhan reservoir (capacity 800 million m³, impounded since 1967) lies in a tectonically active zone of design intensity degree 7 on the MSK-64 scale, and the national monitoring program explicitly tracks changes in the local seismic regime during impoundment and operation. The November 3, 2025 event, however, was a far-field tectonic earthquake originating in Afghanistan ($\mathrm{\Delta }\approx$ 146 km), not a reservoir-triggered event: its hypocentral depth (30 km) and epicentral location are inconsistent with shallow RTS, which typically occurs within a few kilometers of the impoundment and at shallow depth. The InSAR record for 2017-2023 is dominated by seasonal deformation, and no accelerating trend indicative of stress accumulation was identified over that period; individual points nevertheless show cumulative trends (Figs. 28-34), so continued point-wise monitoring remains necessary. The controlled operation of the reservoir – gradual seasonal drawdown and refilling rather than rapid uncontrolled release – further limits any RTS potential. Continuous monitoring nevertheless remains warranted to detect any future change in the local seismic regime.

4.6. Practical implications and comparison with prior studies

The results extend and complement previous studies of earthfill-dam seismic response. Whereas Sultanov and Umarkhonov [2] and Umarkhonov and Loginov [3] analyzed the stress-strain state of earth dams primarily through numerical modeling, and Rashidi et al. [4], Wei et al. [5], and Qi et al. [6] focused on construction-stage or purely numerical dynamic and liquefaction analyses, the present work links such modeling directly to a dual observational dataset (instrumental accelerograms and multi-year InSAR) for a specific recorded event. This observational anchoring is the principal practical advantage: the model predictions are compared quantity-by-quantity with the available field measurements (Table 3), which support the intensity prediction, show the SDOF spectral bound and the measured crest displacement to be not directly comparable quantities, and reveal the PGA gap discussed above.

From an operational standpoint, the demonstrated workflow – automatic earthquake detection on the national electronic platform, retrieval of instrumental accelerograms, decomposition of open Sentinel-1 InSAR data, and rapid evaluation of a transparent chain of analytical models – can be applied to the rapid post-earthquake screening of the 25 reservoir dams affected by the November 3, 2025 event and, more generally, to any monitored dam in a seismically active region. Because the method relies on freely available satellite data and standard analytical models, it offers a low-cost, replicable safety-assessment tool that does not require expensive dedicated instrumentation beyond the existing seismic stations.