Effects of Topographic Factors and Human–Land Relationships on Land‐Use Patterns in the Zhaotong Section of the Jinsha River Basin

1. Introduction

The Xiahaizi landslide in Xiaolongdong, Zhaoyang District, Zhaotong City, is located within the tectonic denudation middle mountain geomorphic region of northeastern Yunnan. The stratigraphic structure is characterized by an unconformable contact between the overlying Quaternary residual slope deposits and the underlying Permian Emeishan basalt formation. The geological conditions are highly complex and structural discontinuities are well developed. Field investigations revealed multiple longitudinal tensile cracks of varying scales distributed along the rear edge of the landslide. These cracks extend for several tens of meters, and their aperture exhibits significant seasonal fluctuations controlled by rainfall and temperature variations. Under the influence of the monsoon climate, the alternating effects of rainfall infiltration during the wet season and freeze thaw cycles during the dry season gradually transform the fracture zone into a preferential groundwater migration pathway, forming a typical seepage stress coupling environment. Localized water accumulation at crack tips increases pore water pressure, inducing hydraulic fracturing and promoting downward crack propagation. Ultimately, the fractures may connect with potential sliding surfaces and become a critical driving mechanism controlling progressive slope failure and overall instability [1,2,3]. Owing to the pronounced heterogeneity and temporal variability of the fracture system during its evolution, conventional borehole investigations and single point monitoring methods are unable to acquire continuous spatial information or effectively characterize the three dimensional structure of crack propagation and damage evolution.

Geophysical methods provide important non destructive technical support for investigating landslide fracture zone structures. Among these methods, the transient electromagnetic method demonstrates significant advantages in identifying water bearing fractures because of its high sensitivity to low resistivity anomalies [4,5]. Previous studies have shown that the transient electromagnetic method possesses strong resolution capability for detecting low resistivity anomalies within collapse zones and landslide bodies, and its response to water saturated fracture zones is superior to that of several conventional electrical methods [6]. Small loop transient electromagnetic systems exhibit high lateral resolution for shallow fracture detection [7], while the apparent magnetic flux method has also been verified as effective for identifying the spatial location of sliding surfaces [8]. However, conventional single phase observation approaches cannot effectively capture the dynamic evolution of electrical parameters during fracture development. In recent years, the development of time lapse geophysical methods has provided new approaches for characterizing time varying subsurface processes. Through comparative analysis of multi temporal datasets, electrical property variations caused by water migration and structural damage can be effectively identified [9,10]. In the field of transient electromagnetics, multiscale inversion methods based on time lapse responses can significantly improve the resolution and stability of resistivity imaging [11]. International studies have systematically investigated the applicability of transient electromagnetic methods in complex terrain regions, fracture identification in volcanic slopes, and anisotropic electrical modeling [12]. Meanwhile, important progress has also been achieved domestically in three dimensional inversion, apparent resistivity algorithms, and multisource information integration [13].

The hydraulic fracturing process within fracture zones is essentially a multiphysical coupling problem involving seepage fields, stress fields, and temperature fields. Existing studies have demonstrated that significant nonlinear coupling relationships exist between seepage flow and normal stress within single fracture rock masses, thereby providing a theoretical basis for understanding fracture propagation [14]. Under water pressure driven conditions, fluids within fractures progressively replace air and form highly conductive pathways, resulting in stage dependent decreases in resistivity [15]. Furthermore, under freeze thaw cycling conditions, water ice phase transitions further intensify fluctuations in electrical parameters, and the corresponding response mechanisms differ substantially from those induced solely by rainfall infiltration [16].

Therefore, constructing a resistivity model for fractured media from an anisotropic perspective and revealing the electrical precursory characteristics of water pressure driven fracture propagation represent key scientific challenges for the quantitative characterization of landslide crack damage. Based on this background, the present study focuses on the Xiahaizi landslide fracture zone in Xiaolongdong, Zhaotong. Time lapse transient electromagnetic high density observations and three dimensional resistivity imaging techniques are employed in combination with field experiments and multitemporal data analysis to systematically investigate the electrical response characteristics during hydraulic fracturing processes. Furthermore, a quantitative relationship between resistivity parameters and fracture damage variables is established to provide methodological support for landslide hazard monitoring and early warning.

2. Geological Background and Fracture Characteristics

The Xiaolongdong Xiahaizi landslide is located in Xiaolongdong Hui and Yi Township of Zhaoyang District in Zhaotong City, Yunnan Province, China. The study area is situated within the transitional zone between the tectonic erosional middle mountain landform region of northeastern Yunnan and the marginal area of the Zhaolu Basin. The region is strongly influenced by the Xiaojiang fault zone and its secondary structures, resulting in complex geological conditions and highly fractured rock masses. The overall topography of the landslide area is characterized by a northwest high and southeast low distribution pattern. The slope gradient is mainly distributed between 15° and 35°, with an average slope angle of approximately 27°. Local steep scarps exhibit slope angles exceeding 40°. The slope surface displays a typical stepped geomorphological structure composed of three to five gentle platforms. The width of each platform ranges from approximately 20 m to 80 m, while adjacent platforms are connected by steep slopes with heights of approximately 5 m to 15 m, forming distinct slope break zones. The elevation of the rear edge of the landslide is approximately 2665 m, whereas the front edge extends to approximately 2500 m, resulting in a relative elevation difference of approximately 165 m. The dominant sliding direction is toward the southeast. A well developed gully system is distributed along the principal sliding direction and exhibits a dendritic pattern. The incision depth of the gullies generally ranges from 5 m to 20 m, with local depths exceeding 25 m. The gully density is approximately 2.0 km/km² to 3.5 km/km². These gullies not only serve as the primary discharge pathways for surface runoff, but also enhance the development and connectivity of internal fractures within the slope mass, thereby providing preferential channels for groundwater seepage and hydraulic fracturing within fracture zones. Previous studies have demonstrated that under similar geomorphological and lithological conditions, the resistivity of fracture zones commonly decreases to 10 Ω·m to 50 Ω·m, indicating the significant development of water bearing fractures.

Figure 1. Geological background and distribution of rear edge tension cracks in the Xiaolongdong Xiahaizi landslide hazard area: (a) Geographic location map of the Xiaolongdong area; (b) Displacement variation derived from InSAR satellite remote sensing monitoring in the Xiaolongdong area; (c) Geomorphological characteristics of the rear edge of the Xiaolongdong landslide; (d) Tension cracks distributed along the rear edge of the Xiaolongdong landslide.

Figure 1. Geological background and distribution of rear edge tension cracks in the Xiaolongdong Xiahaizi landslide hazard area: (a) Geographic location map of the Xiaolongdong area; (b) Displacement variation derived from InSAR satellite remote sensing monitoring in the Xiaolongdong area; (c) Geomorphological characteristics of the rear edge of the Xiaolongdong landslide; (d) Tension cracks distributed along the rear edge of the Xiaolongdong landslide.

The lithological assemblage of the landslide strata provides the material basis for fracture development and hydraulic fracturing processes. The exposed strata within the study area mainly consist of Quaternary residual colluvial deposits and the underlying Permian Emeishan basalt formation. The Quaternary residual colluvial layer is widely distributed from the middle section to the frontal part of the slope, with a typical thickness ranging from 5 to 25 m and locally exceeding 30 m, representing a typical landslide accumulation structure. This layer is primarily composed of brownish red to yellowish brown silty clay, with a clay content generally greater than 30%. The natural water content commonly ranges from 20% to 35% and may exceed 40% under locally saturated conditions. Approximately 15% to 30% basalt breccia fragments are embedded within the soil matrix. These fragments are mostly angular, poorly sorted, and range in size from 2 to 50 mm. Such a structure exhibits pronounced heterogeneity and distinct structural weak surfaces. The permeability coefficient is generally within the range of 10⁻⁶ to 10⁻⁵ m/s. Significant softening occurs after water infiltration, and the shear strength parameter represented by cohesion may decrease by more than 30%. In addition, the relatively high content of clay minerals, including montmorillonite and illite, increases the sensitivity of the material to variations in water content, thereby promoting swelling shrinkage deformation and fracture propagation.

The spatial distribution of landslide fractures is jointly controlled by regional tectonic joints and geomorphological conditions, exhibiting clear zonal development characteristics. A major tensile fracture is developed at the rear edge of the slope, and its orientation is generally consistent with the steep scarp at the slope crest. The overall strike direction is approximately NE50°, and the continuous surface exposure length is about 120 m. Field investigations and analogous landslide cases indicate that the aperture of the rear edge fracture commonly varies significantly, generally ranging from 0.3 to 2.0 m and locally exceeding 3 m. In the study area, the middle section of the fracture exhibits a relatively larger opening width, whereas both ends gradually taper out, which is consistent with the typical distribution pattern of traction induced landslide fractures characterized by a wide central section and narrow margins. The visible depth of the fracture is generally between 3 and 5 m. In some locations, dissolution processes and gravitational collapse have resulted in bead shaped collapse pits or sinkhole like structures. These characteristics are consistent with the coupled development pattern of landslide fractures and seepage observed in southwestern China.

Figure 2. Stratigraphic lithology and tensile cracks in the landslide area: (a) Field investigation of rear edge tensile cracks; (b) Stratigraphic lithology of the landslide area.

Figure 2. Stratigraphic lithology and tensile cracks in the landslide area: (a) Field investigation of rear edge tensile cracks; (b) Stratigraphic lithology of the landslide area.

The middle section of the slope is characterized by the development of multiple groups of secondary fractures. Most of these fractures intersect the main fractures at small angles or exhibit feather like distributions. Their extension lengths are generally between 20 and 50 m, with apertures ranging from 0.1 to 0.3 m. The fractures mainly develop along pre existing structural joint surfaces and commonly display zigzag or en echelon distribution patterns. Transverse tensile fractures dominate the frontal part of the slope. Their orientations are nearly perpendicular to the slope direction, with lengths generally ranging from 10 to 30 m and relatively small apertures. These characteristics indicate deformation controlled by the combined effects of shear compression and tensile stress at the slope front. In cross sectional profiles, the extension depth of fractures is strongly controlled by lithological interfaces and weathering zonation. Most fractures gradually close or transform into bedding parallel shear structures after penetrating into strongly weathered bedrock layers. This phenomenon is consistent with the propagation inhibition effect of fractures in heterogeneous media.

Figure 3. Rain fall statistics chart of Xiaolong Cave.

Figure 3. Rain fall statistics chart of Xiaolong Cave.

The groundwater system within the landslide area is primarily controlled by atmospheric precipitation infiltration recharge, and its hydrological response exhibits significant seasonal and event driven characteristics. The study area is characterized by a typical subtropical monsoon climate with an average annual precipitation of approximately 740 mm. More than 65% of the annual rainfall occurs between June and September, during which intense rainfall events are highly concentrated and of short duration. The rainfall intensity during a single storm event may reach 30 to 80 mm/d. Previous studies have demonstrated that under similar climatic conditions, rainfall infiltration can trigger a rapid increase in slope pore water pressure within several hours, while the lag time of groundwater system response to precipitation is generally limited to several hours or less than one day. During the infiltration process, part of the precipitation rapidly penetrates through fractures and macropores to form deep groundwater recharge, whereas the remaining portion is converted into surface runoff and discharged into the gully system, resulting in a coupled infiltration runoff recharge mechanism. Within the residual slope deposit, the pore permeability coefficient generally ranges from 10⁻⁶ to 10⁻⁵ m/s, leading to relatively slow groundwater migration. In contrast, fractures and structural joints constitute highly permeable pathways with equivalent permeability coefficients ranging from 10⁻⁴ to 10⁻³ m/s, allowing groundwater to respond significantly within several hours after rainfall events. Overall, groundwater flow is directed downslope from the rear edge toward the front edge of the slope. However, within fracture intensive zones, the groundwater flow field undergoes substantial reorganization and exhibits a pronounced convergence toward fracture zones.

Figure 4. Temperature and humidity variations within the tensile crack zone at the rear edge of the Xiaolongdong landslide.

Figure 4. Temperature and humidity variations within the tensile crack zone at the rear edge of the Xiaolongdong landslide.

Zhaoyang District is located in a high altitude mountainous region in northeastern Yunnan Province, where winter temperatures frequently decrease below 0 ℃ and freeze thaw processes are highly active. Meteorological records indicate that from October to December the number of days with mean daily temperatures fluctuating around 0 ℃ exceeds 20 d annually, while the duration of individual freeze thaw cycles generally ranges from 2 to 7 d. These conditions represent a typical seasonal freeze thaw environment. Previous studies have demonstrated that pore water undergoes volumetric expansion of approximately 9% during freezing, generating significant frost heaving pressure within confined fractures. Experimental investigations and numerical simulations indicate that ice wedge pressure within open fractures generally ranges from 5 to 10 MPa, whereas under semi confined or confined conditions the pressure may increase to 10 to 60 MPa and can exceed 100 MPa under extreme conditions. Such stress levels approach or exceed the tensile strength of most weathered rock masses and are therefore capable of promoting crack tip propagation.

Extensive fracture systems are developed within the landslide masses in Zhaoyang District. Residual moisture retained within fractures after the rainy season undergoes repeated freezing during winter, producing persistent frost wedging effects. At the macroscopic scale, monitoring results indicate that fracture apertures exhibit periodic fluctuations during freeze thaw seasons. Compared with autumn conditions, fracture opening in winter generally increases by approximately 5% to 15%. Although partial shrinkage occurs after spring thawing, irreversible plastic deformation is commonly observed. Under the cumulative influence of multiple freeze thaw cycles, major fractures exhibit continuous propagation. In typical landslides, the trailing edge fractures advanced approximately 6 to 10 m within a five year period.

At the microscopic scale, freeze thaw cycling causes significant deterioration of rock and soil structures. Experimental studies have shown that the cohesion of cohesive soils decreases by approximately 20% to 40% after repeated freeze thaw processes, accompanied by substantial increases in porosity and microcrack density. In hard rocks such as basalt, freeze thaw action induces fatigue damage along pre existing structural planes, thereby enhancing fracture connectivity and reducing overall mechanical strength. In addition, freeze thaw processes significantly increase fracture permeability, creating favorable pathways for subsequent rainfall infiltration. During the spring thawing period, the enlarged fracture spaces and weakened fracture walls become more susceptible to hydraulic fracturing under water pressure.

3. Materials and Methods

The transient electromagnetic data acquisition was conducted using an integrated dual coil 600D system in which the transmitting and receiving coils were arranged in a coincident configuration with a coil diameter of 1.2 m. The system operated under a zero offset configuration, which effectively eliminated the blind zone caused by signal loss in the early time channels commonly observed in conventional large loop systems within the time range of 10⁻⁵ to 10⁻⁴ s. This configuration significantly improved the detectability of shallow subsurface responses. Under identical current conditions, the small loop system with a characteristic scale of approximately 1 m exhibited higher resolution for low resistivity anomalies within the depth range of 0 to 30 m. Its equivalent investigation depth and temporal window coverage could extend from several tens of meters to nearly one hundred meters. The dual coil configuration, based on the coaxial arrangement of transmitting and receiving coils, enhanced the amplitude of early induced electromotive force signals and improved the signal to noise ratio. This system is particularly suitable for the high resolution investigation of shallow low resistivity water bearing structures such as landslide fracture zones.

The transmitted current waveform is a bipolar trapezoidal wave with a fundamental frequency of 25 Hz. The turn off time varies with the magnitude of the transmitted current. When the transmitted current is 12 A, the measured turn off time is approximately 2.8 μs. Accurate measurement of the turn off time is critical for late time signal correction. Before each measurement, a no load calibration of the turn off time is performed to ensure the accuracy of the excitation waveform representation in subsequent inversion. The transmitted magnetic moment is calculated as follows:

M = I × N × A

(1)

where I denotes the transmitted current, N represents the number of coil turns, and A is the effective coil area. Substituting I equal to 12 A, N equal to 6 turns, and A equal to 1.13 m² yields a transmitted magnetic moment of 81.4 A m².

The receiver system employs 24 bit analog to digital conversion with a sampling rate of 1 MHz. The sampling window ranges from 2 μs to 8 ms, covering the full time span from early to late stage signals. To address 50 Hz power frequency interference in the landslide area, a notch filter is applied at the hardware level to attenuate the power frequency component. At the software level, signal stacking is performed using an integer multiple of a 20 ms period, with 64 stacks implemented to effectively suppress random noise.

Transient electromagnetic signals exhibit a wide dynamic range, where early time amplitudes reach the volt level and late time amplitudes decrease to the microvolt level. The receiver system adopts a programmable gain amplifier that automatically adjusts amplification according to the sampling stage, with a gain increase of 40 dB for late time signals. Before each measurement, a system self check is conducted to monitor coil resistance and capacitance parameters, ensuring consistency of the coil state.

Figure 5. These figures show the dual coil 600D integrated transient electromagnetic system: (a) 600D transient electromagnetic system equipment; (b) Schematic diagram illustrating the principle of transient electromagnetic detection of anomalous bodies.

Figure 5. These figures show the dual coil 600D integrated transient electromagnetic system: (a) 600D transient electromagnetic system equipment; (b) Schematic diagram illustrating the principle of transient electromagnetic detection of anomalous bodies.

Distortion of the induced electromotive force during the turn off period is corrected using a ramp correction method. The measured turn off waveform is convolved with the inverse Green function to eliminate the influence of the turn off effect on early time data

V_c(t)=V_m(t) × g⁻¹(t)

(2)

where V_c(t) denotes the corrected response, V_m(t) denotes the measured response, and g (t) represents the turn off waveform function.

The Xiaolongdong Xiahaizi landslide area is characterized by significant topographic variation and moderate vegetation coverage. The survey layout must balance terrain adaptability with spatial sampling density. The survey area extends from the main tensile crack zone at the rear edge of the landslide to the shear outlet at the front edge. Ten survey lines are arranged along the slope direction with a spacing of 1 m. Three transverse connecting survey lines are deployed perpendicular to the slope direction, forming a grid observation system. The point spacing is refined to 0.5 m in densely fractured zones and relaxed to 5 m in non fractured areas. A total of 187 physical measurement points are deployed.

Figure 6. Layout of transient electromagnetic survey points and inversion grid constrained by topography.

Figure 6. Layout of transient electromagnetic survey points and inversion grid constrained by topography.

Each measurement point is surveyed using a real time kinematic differential global positioning system to obtain three dimensional coordinates. The horizontal positioning error is less than 0.05 m and the elevation error is less than 0.1 m, providing accurate surface geometric information for subsequent terrain constrained inversion. In areas with complex topography, coil orientation control significantly affects data quality. At each measurement point, a circular area with a diameter of 2 m is manually leveled, and shrubs and loose soil are removed to ensure that the coil plane remains horizontal. For sites with slopes greater than 20 degrees, a tripod is used to support the coil, and dual axis bubble levels are applied to adjust the coil inclination angle to less than 3 degrees. Coil orientation parameters including azimuth and dip angle are recorded in the field for each measurement point and used for orientation correction in subsequent response calculations. When the coil is tilted, the effective receiving area decreases and the induced electromotive force is proportional to the cosine of the tilt angle

V(θ)=V₀ × cosθ

(3)

where V₀ is the induced electromotive force under horizontal conditions and θ is the angle between the coil plane and the horizontal plane. Data acquisition timing is synchronized with environmental conditions. Ambient temperature humidity and wind speed at each measurement time are recorded in a metadata file for reference during data preprocessing. Three repeated observation points are established around the survey area. After every 20 measurement points, repeated observations are conducted to monitor system drift and the influence of external environmental variations on the data. Statistical analysis of repeated observations shows that the relative deviation of induced electromotive force between repeated measurements at the same location is less than 3 percent, indicating that system stability meets the required standards.

Three dimensional inversion is implemented using an unstructured tetrahedral mesh to accommodate complex terrain and irregular anomalous bodies. The terrain surface is constructed based on measured elevation data using triangular meshing and extended downward to generate tetrahedral elements. Near surface elements are refined to smaller sizes, with resolution increased to 0.5 m in fracture zones, while deeper elements are gradually coarsened to 5 m. The total number of elements is approximately 150000, which ensures near surface resolution while controlling computational cost. Forward modeling is based on the time domain finite element method to solve the electric field diffusion equation

$$\nabla \times \left({\displaystyle \frac{1}{\mathsf{\mu }}}\nabla \times \mathrm{E}\right)+\mathsf{\sigma }{\displaystyle \frac{\text{∂}\mathrm{E}}{\text{∂}\mathrm{t}}}=-{\displaystyle \frac{{\text{∂}\mathrm{J}}_{\mathrm{s}}}{\text{∂}\mathrm{t}}}$$

(4)

where E is the electric field intensity, μ is magnetic permeability, σ is electrical conductivity and J_s is the source current density. Time discretization adopts a backward Euler scheme and the solution is advanced from early time to late time steps. The inversion objective function consists of a data misfit term and model constraint terms

$$\mathsf{\Phi }\left(\mathrm{m}\right)=\Vert \mathrm{d}-\mathrm{F}\left(\mathrm{m}\right){\text{‖}}_2^2+\mathsf{\lambda }\parallel {\mathrm{W}}_{\mathrm{m}}\mathrm{m}{\text{‖}}_2^2+\mathsf{\eta }{\Vert {\mathrm{W}}_{\mathrm{s}}\mathrm{m}\Vert }_1$$

(5)

where m is the model parameter vector, d is the observed data vector, F is the forward operator, λ is the regularization parameter, Wm is the smoothness constraint matrix, η is the sparsity weighting factor and Ws is the wavelet transform matrix. The smoothness constraint uses the L2 norm to suppress spurious oscillations in the inversion results and to ensure a reasonable lateral gradient in resistivity distribution. The wavelet sparsity constraint uses the L1 norm to enhance local anomalies such as fracture zones and to avoid boundary blurring caused by over smoothing. The regularization parameter follows a cooling strategy, starting with a large value to ensure rapid convergence toward the background model and gradually decreasing to reveal fine scale features. The inversion is solved using the Gauss Newton iterative method. In each iteration the sensitivity matrix is computed and the regularized linear system is solved. The iteration stops when the data misfit change is less than 1 percent or when the maximum number of iterations reaches 20.

4. Results and Analysis

The field monitoring experiment was conducted in the Xiaolongdong Xiahai landslide fracture zone in Zhaoyang District, Zhaotong City. The monitoring system consisted of a time lapse transient electromagnetic survey network, GNSS displacement monitoring stations, an automatic rain gauge and temperature and humidity sensors. The transient electromagnetic survey network included 187 physical measurement points and completed five data acquisition campaigns over a twelve month period. The sampling periods corresponded to the rainy season, dry season and freeze thaw period. Three GNSS monitoring stations were installed on both sides of the main fracture zone with a sampling frequency of 1 Hz. Daily displacement was obtained after data processing. The rain gauge recorded hourly precipitation, while the temperature and humidity sensors recorded air temperature and relative humidity. During the monitoring period, the total precipitation was 840.37 mm. The precipitation from August to December was 573.85 mm, accounting for 68.2 percent of the total rainfall during the study period. A typical heavy rainfall event occurred from 15 to 18 August with a cumulative precipitation of 86 mm and a maximum hourly rainfall intensity of 24 mm per hour.

Figure 7. Time series diagram of rainfall variation.

Figure 7. Time series diagram of rainfall variation.

Rainfall infiltration is the primary driving force for hydraulic fracturing in fractures. Precipitation preferentially infiltrates rapidly along fracture pathways and accumulates at fracture tips, leading to the development of pore water pressure. When the pore water pressure exceeds the critical threshold for fracture propagation, crack initiation and extension occur. The inversion results shown in Fig. 9 indicate that the electrical resistivity of the fracture zone is highly sensitive to rainfall events. Taking a representative profile of the main fracture as an example, during the early dry season the resistivity ranges from 180 to 250 Ω·m, and low resistivity anomalies appear as isolated clusters. At the end of the dry season, resistivity slightly increases to 200 to 280 Ω·m, reflecting residual ice effects induced by winter freeze thaw processes. In the middle of the rainy season, resistivity sharply decreases to 45 to 85 Ω·m, with a reduction of 55% to 70%, and the low resistivity zones expand and become interconnected along the fracture trend, forming a continuous belt approximately 15 m in width. At the end of the rainy season, resistivity slightly rebounds to 60 to 100 Ω·m but remains significantly lower than dry season levels. At the beginning of the dry season in the following year, resistivity recovers to 120 to 180 Ω·m, indicating partial self sealing effects, although values remain 20% to 40% lower than those of the same period in the previous year, suggesting cumulative irreversible damage.

The influence of freeze thaw cycles on electrical response is mainly reflected during the winter low temperature period. Prior to data acquisition, the study area experienced a freeze thaw season from February to March. During this period, the resistivity of the fracture zones increased by 15 percent to 25 percent compared with the previous stage, while GNSS displacement data indicated a slight increase in fracture aperture of approximately 0.5 mm. This apparent inconsistency can be explained by the high resistivity characteristics of ice. When pore water within fractures freezes into ice, the resistivity increases sharply from 10 to 30 ohm meter for liquid water to values exceeding 10 to the power of 4 ohm meter, resulting in an overall increase in apparent resistivity. After ice melting in spring, meltwater temporarily recharges the fractures, leading to a rapid decrease in resistivity. The post thaw resistivity values are lower than those observed before freezing, indicating irreversible fracture expansion induced by frost heave.

Figure 8. Comparison of resistivity inversion results of fracture zones at different time periods: (a)(b)2025-5-11;(c)2025-5-19;(d)®2025-6-4;(f)(g)2025-6-6;(h)2025-6-15;(i)2025-6-16;(j)2025-6-23;(k)2025-6-26;(l)2025-6-28;(m)(n)2025-7-26;(o)(p)2025-7-29;(q)®2025-8-5.

Figure 8. Comparison of resistivity inversion results of fracture zones at different time periods: (a)(b)2025-5-11;(c)2025-5-19;(d)®2025-6-4;(f)(g)2025-6-6;(h)2025-6-15;(i)2025-6-16;(j)2025-6-23;(k)2025-6-26;(l)2025-6-28;(m)(n)2025-7-26;(o)(p)2025-7-29;(q)®2025-8-5.

Coupled analysis shows that regions undergoing freeze thaw cycles exhibit a 10 percent to 15 percent greater reduction in resistivity during the rainy season compared with regions not affected by freeze thaw processes. Approximately 30 percent of the annual cumulative displacement is attributed to direct contributions from the freeze thaw period. Field monitoring experiments further reveal the complete electrical response process of rainfall induced hydraulic fracturing in fractures, characterized by an initial slight increase in resistivity, followed by accelerated decline, low resistivity connectivity, and recovery during pressure release. A quantitative relationship is established between the rate of resistivity change and fracture aperture. These results provide a physical basis and interpretative framework for monitoring damage evolution in landslide fractures.

5. Discussion

The GNSS displacement data and TEM inverted resistivity obtained from field monitoring tests provide a basis for establishing a quantitative relationship. Three GNSS monitoring points on both sides of the main fracture were selected. The relative displacement between points on opposite sides of the fracture was calculated at each data acquisition time to estimate the average fracture aperture d. The corresponding TEM inverted resistivity values at the same locations were extracted to construct a dataset consisting of five monitoring periods and forty samples. The resistivity was normalized, and the resistivity variation rate was defined as

Δρ = (ρ_₀-ρ) / ρ_₀

(6)

where ρ₀ represents the background resistivity measured in the first period during the initial dry season and ρ represents the measured resistivity at a given time.

Based on the parallel conduction model, the equivalent resistivity of fractured rock mass with water filled fractures exhibits a theoretical relationship with fracture aperture. The derived approximate expression can be reformulated as

$$\mathsf{\Delta }\mathsf{\rho }=\frac{{\mathsf{\rho }}_{\mathtt{r}}·\mathsf{\rho }}{{\mathsf{\rho }}_{\mathtt{r}}}=·\frac{\frac{{\mathsf{\rho }}_{\mathsf{\omega }}}{{\mathsf{\rho }}_{\mathtt{r}}}\frac{\mathtt{d}}{\mathrm{L}}}{1-\frac{{\mathsf{\rho }}_{\mathtt{r}}}{{\mathsf{\rho }}_{\mathsf{\omega }}}\frac{\mathtt{d}}{\mathrm{L}}}.$$

(7)

where ρ_r is the resistivity of the rock matrix, ρ_w is the resistivity of fracture water, d is the fracture aperture and L is the representative elemental length. This expression has a rational function form. When d is small, Δρ is approximately linear with respect to d. When d increases, Δρ approaches unity.

Field data fitting was performed using nonlinear least squares, with the objective function defined as the sum of squared residuals between predicted and observed values. The results indicate that the parameter ρr over ρw is 185 and L is 3.2 m. The coefficient of determination R squared is 0.92 and the root mean square error is 0.038. The fitted curve agrees well with observed data when the aperture ranges from 0.2 mm to 2.8 mm. When the aperture is less than 0.2 mm, the fitted values are slightly higher than observations, indicating a delayed resistivity response during the micro fracture stage. When the aperture exceeds 2.8 mm, the observed values are slightly lower than the fitted curve, which may be attributed to reduced contact points on fracture walls at high aperture conditions, leading to more efficient conductive pathways than those assumed in the idealized model.

Figure 9. Field measured resistivity data points and fitted curves of tensile fractures at different stages in the rear boundary of the Xiaolongdong landslide.

Figure 9. Field measured resistivity data points and fitted curves of tensile fractures at different stages in the rear boundary of the Xiaolongdong landslide.

The data distribution indicates that the resistivity variation rate exhibits a nonlinear increase with increasing fracture aperture. When the aperture is less than 0.5 mm, the growth of Δρ is relatively slow. As the aperture increases from 0 to 0.5 mm, Δρ ranges approximately from 0.12 to 0.20. In the aperture interval of 0.5 to 2.5 mm, the growth rate of Δρ accelerates, with each 0.5 mm increase in aperture corresponding to an increase in Δρ of approximately 0.15 to 0.25. When the aperture exceeds 2.5 mm, Δρ gradually approaches saturation, with a maximum value of approximately 0.75 to 0.85, and further increases in aperture produce limited changes in resistivity. This trend reflects the process of conductive pathway formation after water saturation in fractures. At the early stage with small apertures, poor hydraulic connectivity leads to limited reduction in resistivity. At the intermediate stage, increasing aperture promotes the formation of continuous high conductivity pathways, resulting in a rapid decrease in resistivity. At the late stage, when the aperture becomes sufficiently large, current flow is mainly governed by the water phase, while the contribution of the surrounding rock becomes negligible, leading to a stabilized resistivity response.

The fracture zone of the Xiaolongdong Xiahaizi landslide develops multiple dominant orientations, resulting in anisotropic characteristics in three dimensional resistivity. The resistivity obtained from transient electromagnetic inversion represents an equivalent isotropic property and cannot accurately describe the control effect of fracture orientation on current flow. Therefore, an anisotropic resistivity model considering fracture development direction is constructed, where resistivity is expressed in tensor form as a second order tensor.

Figure 10. GNSS displacement data validation of the tensile crack evolution at the rear edge of the Xiaolongdong landslide after failure: (a) cumulative horizontal displacement field;(b) cumulative vertical displacement field;(c) temporal evolution of crack aperture derived from GNSS monitoring;(d) spatial distribution of TEM-derived resistivity anomalies;(e) correlation between displacement and resistivity variation rate;(f) validation of crack damage evolution inferred from TEM inversion using GNSS observations.

Figure 10. GNSS displacement data validation of the tensile crack evolution at the rear edge of the Xiaolongdong landslide after failure: (a) cumulative horizontal displacement field;(b) cumulative vertical displacement field;(c) temporal evolution of crack aperture derived from GNSS monitoring;(d) spatial distribution of TEM-derived resistivity anomalies;(e) correlation between displacement and resistivity variation rate;(f) validation of crack damage evolution inferred from TEM inversion using GNSS observations.

The calibrated model is used to compute the principal axis resistivity and the anisotropy coefficient at different positions within the fracture zone. The principal axis resistivity is defined as the reciprocal of the three eigenvalues of the conductivity tensor, corresponding respectively to the directions of maximum, intermediate and minimum resistivity. The anisotropy coefficient λ is defined as the ratio between the maximum and minimum resistivity. The results indicate that λ within the fracture zone ranges from 2.5 to 6.8. Higher values of λ are observed in the densely fractured rear region, reaching a maximum of 6.8, whereas lower values occur in the sparsely fractured frontal region with a minimum of approximately 2.5. The direction of maximum resistivity is approximately perpendicular to the main fracture planes, while the direction of minimum resistivity is parallel to the fracture strike, which is consistent with theoretical expectations.

The anisotropic model modifies the resistivity values obtained from isotropic inversion. The isotropic results fall between the maximum and minimum resistivity values and tend to represent a weighted average, which fails to accurately capture the true current flow paths. After applying the anisotropic model, the spatial correlation between inversion results and GNSS displacement increases from 0.61 to 0.78, indicating improved consistency with real conditions.

Figure 11. Relationship between damage variable D and normalized apparent resistivity index Rn.

Figure 11. Relationship between damage variable D and normalized apparent resistivity index Rn.

Table 1. Statistical characteristics of damage variable D for fracture zones in different stages.

Monitoring period	Range of D values in rear main fracture zone	Range of D values in central secondary fracture zone	Range of D values in frontal fracture zone	Range of D values in frontal fracture zone Corresponding relationship with GNSS displacement
Early dry season	<0.20	<0.18	<0.15	Low displacement rate zone
Mid rainy season	0.45-0.58	0.32-0.42	0.25-0.35	Areas with displacement rates greater than 15 millimeters per month show D values above 0.40
Late rainy season	0.35-0.45	0.25-0.35	0.18-0.25	Displacement rate decreases but remains higher than dry season levels

Table 1. Statistical characteristics of damage variable D for fracture zones in different stages.

Monitoring period	Range of D values in rear main fracture zone	Range of D values in central secondary fracture zone	Range of D values in frontal fracture zone	Range of D values in frontal fracture zone Corresponding relationship with GNSS displacement
Early dry season	<0.20	<0.18	<0.15	Low displacement rate zone
Mid rainy season	0.45-0.58	0.32-0.42	0.25-0.35	Areas with displacement rates greater than 15 millimeters per month show D values above 0.40
Late rainy season	0.35-0.45	0.25-0.35	0.18-0.25	Displacement rate decreases but remains higher than dry season levels

Based on the three dimensional inversion results of time lapse Transient Electromagnetic Method data acquired in the Xiahaizi landslide area of Xiaolongdong, the spatiotemporal evolution characteristics of the permeability tensor within the landslide fracture zone were analyzed. During the early dry season, the anisotropy of the fracture zone permeability coefficient was relatively weak, and no distinct pattern was observed in the direction of maximum permeability. The maximum permeability coefficient ranged from approximately 2.5 × 10⁻⁶ to 8.6 × 10⁻⁶ m/s, indicating that the fractures remained in a relatively closed state and that hydraulic fracturing effects were limited.

During the middle stage of the rainy season, continuous precipitation increased fracture water pressure and intensified hydraulic fracturing processes, resulting in a significant increase in fracture aperture. The overall permeability coefficient increased by approximately 1.5 orders of magnitude, while the maximum permeability coefficient increased to 1.2 × 10⁻⁴ to 6.8 × 10⁻⁴ m/s. Permeability anisotropy became substantially enhanced. The direction of maximum permeability was consistent with the orientation of the major fractures, whereas the minimum permeability direction was perpendicular to the fracture planes. These observations demonstrate the dominant control exerted by hydraulic fracturing on seepage flow within the fracture network.

During the late rainy season, as precipitation intensity decreased and geostress was partially released, partial closure of fractures occurred. Consequently, the permeability coefficient declined to 3.5×10⁻⁵ to 1.2×10⁻⁴ m/s. However, the maximum permeability coefficient remained higher than that observed during the early dry season, suggesting the existence of irreversible permeability enhancement induced by hydraulic fracturing. The evolution trend of the permeability tensor obtained from time lapse TEM inversion exhibited strong correlations with surface precipitation and groundwater level variations. Following a single intense rainfall event, the peak permeability coefficient appeared with a delay of approximately two to three days, which was consistent with the delayed response of groundwater levels. This phenomenon indicates a pronounced hysteresis effect in the hydrodynamic response of fractures.

Seasonal freeze thaw processes also produced periodic fluctuations in fracture permeability coefficients, although the magnitude of variation was significantly lower than that observed during the rainy season. The maximum permeability coefficient generally ranged from 3×10⁻⁶ to 1.1×10⁻⁵ m/s, indicating that the influence of freeze thaw processes on fracture permeability was secondary compared with hydraulic fracturing effects.

6. Conclusions

In this study, time-lapse three-dimensional transient electromagnetic (TEM) imaging was employed to investigate the seasonal evolution of hydraulic fracturing and damage processes within the Xiaolongdong–Xiahaizi landslide fracture zone. By integrating multi-temporal geophysical observations with deformation monitoring and hydrogeological analyses, a comprehensive framework was established to characterize the spatiotemporal evolution of fracture systems and to evaluate their implications for landslide instability.

The results reveal pronounced seasonal variability in the electrical structure of the landslide fracture zone. During the early dry season, the fracture zone is characterized by relatively high apparent resistivity values ranging from 180 to 320 Ω·m, while low-resistivity anomalies occur as isolated discontinuous patches. Following rainfall infiltration during the wet season, apparent resistivity decreases by approximately 55–70%, and low-resistivity anomalies become progressively interconnected along the northeast direction, with resistivity values declining to 45–95 Ω·m. The most significant electrical anomalies are concentrated at depths of 5–15 m, corresponding closely to the identified sliding surface. During freeze–thaw periods, resistivity increases by 15–25% within 3–8 m of the fracture boundaries, reflecting the influence of frost-heave-induced structural changes. Furthermore, the spatial extent of resistivity variation exhibits a strong positive correlation with precipitation intensity. GNSS observations indicate that cumulative displacement in areas exceeding the identified resistivity threshold is approximately 2.3 times greater than that in surrounding stable zones, demonstrating a strong coupling relationship between electrical property evolution and landslide deformation.

Quantitative analyses further demonstrate that the resistivity variation rate (Δρ) is highly correlated with fracture aperture development, yielding a coefficient of determination of 0.93 and enabling reliable inversion of fracture apertures ranging from 0.1 to 2.5 mm. The resistivity anisotropy ratio varies between 2.3 and 7.1, with the principal resistivity direction consistently oriented perpendicular to the dominant fracture plane. A significant linear relationship is identified between the damage variable (D) and normalized apparent resistivity (Rn), expressed as D = 0.38Rn + 0.09, providing a quantitative indicator of early-stage fracture damage. Permeability tensor inversion further indicates that the maximum permeability coefficient reaches 10⁻⁴–10⁻³ m·s⁻¹ during the wet season, representing an increase of one to two orders of magnitude relative to dry-season conditions. The principal permeability direction coincides with the dominant fracture orientation, while the permeability response lags rainfall events by approximately 2–3 days, highlighting the delayed hydromechanical response of the landslide system.

The time-lapse TEM approach demonstrates strong capability for imaging shallow fracture networks within depths of approximately 20 m, achieving horizontal positioning errors of less than 1 m and vertical errors of less than 2 m. Fractures with apertures exceeding 0.3 mm can be effectively detected and monitored. Distinct precursor signals characterized by an initial resistivity increase followed by a rapid decrease of 3–8% are consistently observed prior to fracture propagation, with inflection points occurring approximately 2–5 min before crack expansion. During unloading stages, partial fracture self-sealing results in residual resistivity reductions of 20–40%. The contrasting electrical responses associated with freeze–thaw processes and rainfall infiltration further demonstrate the sensitivity of TEM observations to different damage mechanisms. Overall, the proposed time-lapse three-dimensional TEM framework provides an effective tool for quantifying fracture evolution, evaluating landslide stability, and supporting early-warning systems in geologically active mountainous regions.