A goaf-slope system stability analysis method based on dynamic strength reduction method

Abstract

A kind of goaf-slope system stability analysis method based on dynamic strength reduction method, comprising the following steps;(a) initial overall strength reduction calculation is carried out to goaf-slope system model, and initial damage area is identified;(b) the mechanical parameters of the starting point area are reduced and iterative calculation is carried out, and the main control disaster source is determined using minimum critical local reduction coefficient criterion;(c) the plastic zone extension path triggered by it is tracked, and the newly born plastic damage zone is gradually included in the dynamic reduction range until the macroscopic failure surface is formed, simulating the complete progressive failure process;(d) introduce zoned progressive dynamic strength reduction strategy, control the attenuation of rock mass strength parameters through field variable, establish the functional relationship between strength reduction coefficient Fv and analysis step time;(e) real-time call multi-parameter fusion instability criterion to determine the limit safety factor of slope.The present application can realize the accurate identification of disaster source in complex rock mass system and the dynamic evaluation of failure process.

Description

Technical Field

[0001] This invention relates to the field of geotechnical engineering stability analysis technology, specifically to a stability analysis method for a goaf-slope system based on the dynamic strength reduction method. Background Technology

[0002] In the process of converting open-pit mines to underground mines or vice versa, complex rock and soil systems often form, where goafs and slopes coexist. The stability of such systems directly affects not only mine safety and resource recovery efficiency but also the safety of the surrounding ecological environment, thus remaining a core concern in mining engineering and geotechnical engineering. Traditional stability analysis methods, such as the limit equilibrium method, are ineffective in handling complex geometries and heterogeneous rock masses. While the global strength reduction method, widely used in numerical simulations, can provide an overall safety factor for the system, its basic idea is to synchronously and proportionally reduce all mechanical parameters of the rock and soil in the overall model. This globally synchronous reduction method has significant limitations: it implicitly assumes that all units within the rock mass system simultaneously reach a limit equilibrium state, while in actual engineering, the failure process is often triggered first by local weak points and gradually expands to the surrounding area under stress redistribution, exhibiting a clear progressive evolution characteristic. Furthermore, traditional strength reduction methods can usually only identify the potential sliding surface that eventually forms, but they struggle to reveal the initiation location, triggering mechanism, and evolution path of failure, lacking an effective characterization of the temporal and causal relationships of the failure process.

[0003] Therefore, when analyzing the chain-like disaster problem of the goaf-slope coupling system, this method often shows that the results are too conservative or fail to reflect the real physical mechanism, thus making it difficult to achieve accurate identification and effective early warning of key disaster-causing areas.

[0004] In summary, existing stability analysis methods are insufficient in revealing the progressive failure mechanism of complex rock mass systems and accurately locating disaster sources, making it difficult to meet the needs of modern mines for refined and scientific safety management. Summary of the Invention

[0005] To overcome the above technical problems, the present invention aims to provide a stability analysis method for goaf-slope systems based on dynamic strength reduction. This method can accurately capture the main disaster-causing source. By constructing a disaster-causing source determination mechanism for goaf-slope systems based on local instability sensitivity analysis, and by comparing the critical local reduction coefficients required for different initial damage zones to trigger overall system instability, the main disaster-causing source can be quantitatively identified. Furthermore, the chain failure mechanism driven by the source can be simulated to clarify the temporal nature of the failure, providing technical support for disaster mechanism analysis and stability evaluation of complex rock mass engineering.

[0006] To achieve the above objectives, the technical solution adopted by the present invention is as follows: A stability analysis method for a goaf-slope system based on dynamic strength reduction includes the following steps; (a) Initial damage zone identification: By performing initial overall strength reduction calculations on the goaf-slope system model, one or more initial damage areas that have reached the plastic yield state within the system are identified. (b) Identification of the main source of disaster: Starting from each of the initial damage zones, local strength reduction is initiated. By iteratively reducing the mechanical parameters of the starting zone and performing calculations, the expansion trend of plastic failure is analyzed. If a certain initial damage zone triggers a chain failure effect and dominates overall instability during the local strength reduction process with the minimum critical local reduction coefficient Fi, then the minimum Fi is used as the basis for determining the main disaster source, and the corresponding initial damage zone is determined as the main disaster source. (c) Progressive Destruction Path Tracing: Taking the main disaster-causing source as the core, the expansion path of the plastic zone caused by it is continuously tracked, and the newly formed plastic failure zone is gradually included in the dynamic reduction range until a continuous macroscopic failure surface is formed, thereby simulating a complete progressive failure process; the expansion process of the dynamic reduction range essentially constitutes a progressive failure propagation mechanism based on spatial adjacency, which is used to characterize the chain-like instability evolution path driven by the main disaster-causing source. (d) Zonal progressive dynamic intensity reduction strategy: As the underlying driving mechanism for parameter reduction in steps (b) and (c), field variables are introduced to control the attenuation of rock mass strength parameters, and a functional relationship between the strength reduction coefficient Fv and the analysis step time is established, so that the cohesion C and internal friction angle φ of the rock mass dynamically decrease as the analysis step progresses until the system reaches the critical instability state. This strategy realizes the non-uniform dynamic evolution of rock mass strength parameters through the coupling of time dimension and spatial partition.

[0007] (e) Criterion for instability of multi-parameter fusion: During the iterative calculation process of steps (b) and (c) above, the multi-parameter fusion instability criterion is called in real time as the loop termination condition. The single non-convergence criterion is abandoned, and a multi-parameter fusion instability criterion based on the sudden change of characteristic point displacement, the plastic zone continuity rate and the sudden increase of system kinetic energy is established. The ultimate safety factor of the slope is determined by monitoring the inflection point of the data curve.

[0008] In step (a), the reduction factor used in the initial overall strength reduction calculation is slightly greater than 1.0, with a value of 1.05-1.2; the criteria for determining the initial damage area are: Under the aforementioned reduction factor, the equivalent plastic strain of the rock mass element reaches a preset threshold, or enters a shear or tensile yield state.

[0009] In step (b), the local strength reduction process is independent and exclusive; for each candidate initial damage zone, after resetting the model state, only the rock mechanics parameters within its region are reduced, including the tangent of cohesion c and internal friction angle φ. The specific steps for calculating the partial reduction are as follows: b1) Set the initial local reduction coefficient Fi = 0.8, the reduction step size ΔF (e.g., 0.05), and the iteration step number i = 0; b2) Calculate the reduction factor for the current iteration step i: F i =0.8+i×ΔF; b3) Based on the current reduction factor F i Calculate the updated mechanical parameters of the rock mass within the initial damaged zone, including the updated cohesion c. i =c / F i The updated internal friction angle φ i =arctan(tanφ / F i ), where c and φ are the initial parameters of the rock mass; b4) Assign the updated parameters to the region cell and run the calculation until equilibrium is reached; b5) Determine whether the system is unstable based on the multi-parameter fusion instability criterion in step (e). If it is not unstable, let i = i + 1 and return to step b2) to continue the loop; if the system is unstable, terminate the loop and record the current F. i As the critical local reduction factor for this region; The initial damage zone with the smallest critical local reduction coefficient is identified as the main source of disaster.

[0010] In step (c), the expansion of the dynamic reduction range follows the trigger-inclusion principle; In one calculation step, intensity reduction is performed only on the determined reduction range (initially the main disaster source); after the calculation is balanced, the units that are adjacent to the current reduction zone and newly enter the plastic state are identified to form a new plastic zone; in the next calculation step, this new plastic zone is merged into the dynamic reduction range, and the next level of intensity reduction is applied to the updated entire reduction range; this process is repeated.

[0011] In step (d), the formula for dynamic intensity reduction is as follows: Among them, C i φ i The reduced cohesion and internal friction angle are represented by C and φ, respectively, which are the initial parameters of the rock mass. Fv is the strength reduction factor. The initial value of Fv is set to 0.8, and it is dynamically increased through the analysis steps.

[0012] In step (e), the criteria for determining the formation of the macroscopic failure surface are: the numerical calculation results do not converge, or the plastic zone extends from the toe of the slope or the boundary of the goaf to the top of the slope or other free surfaces, forming a continuous failure zone; the final stability coefficient of the system is determined by the critical local reduction coefficient of the main disaster source and the reduction in the subsequent progressive failure process.

[0013] In the multi-parameter fusion instability criterion, the weights of each criterion are fused through normalization to improve the stability and reliability of the instability determination.

[0014] The method is implemented on numerical simulation software platforms such as finite element or finite difference (e.g., FLAC3D, ABAQUS). Through embedded FISH language or Python script programming, it realizes dynamic, partitioned, iterative modification of model element parameters and automated control of the calculation process.

[0015] The specific script-based automated control steps are as follows: 1) Model and variable initialization: Establish the geometric model and assign initial mechanical parameters, and define global variables to store the characteristic point displacement, element plastic state and system kinetic energy for each analysis step; 2) Partition identification and parameter traversal module: Write a script to read the plastic strain data of the unit and automatically divide the units that exceed the threshold into independent groups as the initial damage area; 3) Dynamic reduction and adaptive optimization module: Write loop control statements to automatically execute the parameter modification formula of step (b) within the loop body, and automatically reset the model and call the solver; 4) Tracking and Boundary Extension Module: After each calculation substep, the script automatically retrieves the elements adjacent to the current reduction zone. If its State indicator changes to plastic yield, it is dynamically added to the current reduction group via command. 5) Instability detection and output module: Extracts displacement, plastic zone volume and kinetic energy data in real time during the calculation process, calculates their derivative with respect to time (i.e. rate of change), and when the rate of change exceeds the set abrupt change threshold, the script automatically triggers an interrupt command to stop the calculation and outputs the current reduction coefficient as the ultimate safety factor of the slope.

[0016] The method described herein is applicable to the analysis of goaf-slope coupling systems composed of various lithologies and complex structural planes (such as faults and joints), and is particularly suitable for evaluating the long-term stability of goaf areas or the final slope of open-pit mines under steep slopes affected by adjacent goaf activities.

[0017] The beneficial effects of this invention are: (1) Accurate search for disaster source: Compared with the traditional overall reduction method, the present invention identifies and selects multiple candidate areas through the initial damage area identification in step (a), and then performs independent and exclusive local intensity reduction calculations on each candidate area in step (b). The minimum critical local reduction coefficient is used as the criterion to adaptively search for the main disaster source. Through the implementation of the above steps, the method can objectively identify the main disaster source that plays a decisive role in the overall stability from multiple potential instability areas, accurately reveal the location of the disaster source, and transform global fuzzy analysis into local accurate simulation, making the analysis of disaster mechanism more in line with physical reality.

[0018] (2) Simulation and reproduction of the failure process: By tracing the chain failure driven by the main disaster source, this method can dynamically simulate the complete progressive failure process of the rock mass from local yielding to macroscopic penetration, clearly showing the path, sequence and spatiotemporal evolution characteristics of the failure, and realizing the visualization playback of the instability process.

[0019] (3) Achieving precise prevention and control of disasters based on their temporal sequence: Unlike traditional methods that can only provide a single safety factor, this method can clarify the causal chain of disasters. Because step (c) of this invention tracks and records the temporal sequence of the formation of the new plastic zone step by step, the progressive damage process of the damage path is clearly output. Based on the temporal damage path determined in this step, there is no need to blindly reinforce the entire system on a large scale. That is, only local reinforcement is needed for the main disaster-causing source determined in step (b) or the early critical expansion node identified in step (c). This can cut off the damage chain at the bud stage, significantly improving the pertinence and economic benefits of disaster prevention and mitigation measures.

[0020] (4) Wide applicability: The regional dynamic strength reduction strategy proposed in this method solves the problem that the traditional overall reduction method cannot objectively reflect the failure sequence caused by strength differences in multiple structures such as goaf, fault, and slope. The technical logic is universal, based on a general numerical simulation platform, implemented through programming, does not depend on specific software, is easy to promote, and can also be widely applied to the stability evaluation of other complex geotechnical engineering projects with significant spatial differences. Attached Figure Description

[0021] Figure 1 This is a schematic diagram of the analysis process logic of the present invention.

[0022] Figure 2 This is a flowchart illustrating the dynamic strength reduction progressive failure simulation of the present invention.

[0023] Figure 3 This is a schematic diagram of the expansion of the plastic zone in this invention.

[0024] Figure 4 This is a geometric model diagram of the finite element numerical simulation of the goaf-slope system of the present invention.

[0025] Figure 5 This is a stress cloud diagram from the finite element numerical simulation of the present invention.

[0026] Figure 6 This is a diagram showing the evolution of the plastic zone in the finite element numerical simulation of this invention.

[0027] Figure 7 This is a curve showing the displacement-reduction coefficient of the key monitoring point in this invention. Detailed Implementation

[0028] The present invention will now be described in further detail with reference to the accompanying drawings.

[0029] like Figure 1 As shown, this invention discloses a stability analysis method for a goaf-slope system based on the dynamic strength reduction method. The analysis process of this invention is shown in the following schematic diagram. The overall process mainly includes the following parts. Initial damage zone search mechanism: This invention first introduces a low initial overall strength reduction factor Fv to perform a "coarse scan" of the entire system. Physically, this means simulating the regions where the system first responds and enters the nonlinear deformation stage when subjected to slight overload or minor degradation of material parameters. According to the Mohr-Coulomb yield criterion, an element enters the plastic state when the stress state satisfies f(σ) = 0. By setting a plastic strain threshold, these regions that yield earliest, i.e., the initial damage zone, can be precisely delineated in the model. These regions are potential weak points in the system and are the focus of subsequent analysis.

[0030] Principle of identifying the main source of disaster: This step is the core innovation of this invention. Its principle is similar to sensitivity analysis of a system. For the multiple identified initial damage areas (such as A, B, and C), this method performs stress tests on each one.

[0031] The specific operation is as follows: In the numerical model, only reduce the strength parameters (cohesion c, internal friction angle φ) of region A, observe how the system stability responds, and record the local reduction factor Fv required for the system to experience overall instability. A Then, the model was reset, and the same stress test was performed on region B to obtain Fv. B And so on.

[0032] Finally, compare all the obtained local reduction coefficients Fv, and select the coefficient with the smallest value, for example, Fv A The corresponding region A is the primary source of the disaster. This is because it means that even the smallest degree of degradation applied to region A can trigger a domino effect in the entire system, leading to chain-reaction damage.

[0033] Progressive destruction dynamic tracking mechanism: like Figure 2 The diagram shown is a detailed flowchart of the dynamic strength reduction progressive failure simulation of this invention. After identifying the main disaster source, this method enters the "scenario simulation" stage, according to... Figure 2 The illustrated cyclical calculation process starts from the main source of disaster and simulates how its damage expands step by step.

[0034] (1) Spatiotemporal evolution simulation: Using an iterative calculation method, at the k-th step, a reduction coefficient Fv is applied to the currently formed plastic region (initially the main controlling disaster source). k After equilibrium is achieved, the stress redistribution in the rock mass leads to the formation of new plastic units in the adjacent area, resulting in the formation of newly formed plastic zones.

[0035] (2) Dynamic reduction zone update: In step k+1, the newly formed plastic region is merged with the original plastic region to form a larger dynamic reduction range. Then, a reduction coefficient Fv is applied to this updated range. （k+1） (k+1>k).

[0036] (3) Termination Criterion for Damage: Repeating the above process, the dynamic reduction range will continue to expand until a macroscopic damage zone is formed that runs from the weak point inside the system to the free surface, such as... Figure 3 The diagram shown illustrates the expansion of the plastic region in this invention. At this point, the numerical calculations will fail to converge due to the dramatic increase in displacement, indicating that the system has experienced overall instability. The entire process completely reproduces the spatiotemporal evolution path of the failure. Specific implementation examples: This embodiment aims to simulate the impact of the goaf on the stability of the upper steep slope during the transition from open-pit to underground mining using ABAQUS finite element analysis software, and to determine its critical safety factor.

[0038] A. Geometric Model Establishment like Figure 4 The figure shows a simplified model of the goaf-slope system, and the geometric model was simulated using the ABAQUS finite element method. The overall dimensions of the model are set as follows: length 100m, height 78m; slope height 40m, slope angle 45°, width 20m.

[0039] Characteristics of the goaf: Located below the toe of the slope, with a depth (vertical distance from the toe) of 10m, a width of 10m, and a height of 8m. The model was meshed using reduced integral elements, resulting in 9840 elements and 10250 nodes. Mesh refinement was applied around the goaf and in the potential sliding surface area of ​​the slope, with the mesh size controlled at 0.25m × 0.25m.

[0040] B. Material Property Definition Based on the results of deep rock mechanics tests and the Mohr-Coulomb criterion, accurate physical and mechanical parameters were defined for each rock layer in the model. The specific values ​​for the model are as follows: Density ρ = 2450 kg / m³ 3 ; Elastic modulus E = 9.0 GPa; Poisson's ratio ν = 0.25; Cohesion c = 351 kPa; Internal friction angle φ = 28°; C. Boundary conditions and loads; Displacement boundary: Horizontal constraints in the X and Y directions are applied to the left and right sides of the model respectively (U1=0, U2=0), and full constraints in the X, Y, and Z directions are applied to the bottom (U1=U2=U3=0). The top surface and slope of the model are free surfaces.

[0041] Load: Apply a gravitational load, with gravitational acceleration taken as 9.8 N / s². 2 .

[0042] D. Simulation analysis process Step 1: Geostress Balance By using the Geostatic analysis step in ABAQUS, gravity is applied and the initial stress field is balanced to ensure that the model undergoes only minimal elastic deformation before excavation.

[0043] Step 2: Excavation of the goaf: By using the Model Change function to remove the set of cells representing the goaf area, stress redistribution occurs inside the slope after the simulation is completed.

[0044] Results analysis (e.g.) Figure 5 As shown in the figure: a significant tensile stress concentration zone appears on the roof of the goaf, and the tensile stress zone extends towards the slope toe, indicating that mining has had an unloading effect on the slope.

[0045] Step 3: Dynamic intensity reduction: By using predefined field variables, the C and tanφ of all geotechnical materials are gradually reduced from 0.8 times to 2.0 times with the analysis step time (i.e., the reduction factor Fv increases linearly from 0.8 to 2.0).

[0046] E. Results Analysis and Safety Factor Determination: Analysis of displacement mutation characteristics: like Figure 7 As shown in b, this is the curve of the horizontal displacement of the monitoring point at the top of the slope selected in this invention as a function of the reduction coefficient Fv.

[0047] When Fv < 1.1, the displacement curve is flat; When Fv reaches 1.15, the curve shows a clear inflection point, the displacement rate increases sharply, and the tangent slope tends to be vertical. Based on the displacement mutation criterion, the safety factor is determined to be 1.15.

[0048] Analysis of the characteristics of the plastic zone: like Figure 6 As shown, this is a plastic strain contour map at different stages of the reduction process.

[0049] Figure 6 a (Fv=0.9): Local plastic yielding occurs only at the apex of the goaf; Figure 6 b (Fv=1.05): The plastic zone at the toe of the slope expands upward, while the fissures in the goaf extend towards the surface; Figure 6 c (Fv=1.16): A continuous plastic shear zone extends upward from the top of the goaf. Based on the plastic zone penetration criterion, the reduction coefficient corresponding to the critical state is approximately 1.16.

[0050] Overall judgment: Combining the above three criteria, the minimum value among them was taken and mutually verified, and the stability safety factor of the goaf-slope system was finally determined to be 1.15. This result indicates that under the current working conditions of the goaf, the slope is in a basically stable state, but monitoring or local reinforcement of the slope toe and the goaf roof is required.

[0051] Compared to traditional calculations that do not consider the impact of goaf areas, such as Figure 7 As shown in a, the safety factor calculated by the method of the present invention is approximately 1.3, which reduces the slope stability by about 11.5%. This accurately quantifies the reduction effect of underground mining on slope stability and verifies the effectiveness and engineering practical value of the method.

Claims

1. A goaf-slope system stability analysis method based on dynamic strength reduction method, characterized in that, Includes the following steps; (a) By performing initial overall strength reduction calculations on the goaf-slope system model, one or more initial damage areas that have reached the plastic yield state are identified within the system; (b) Starting from each of the initial damage zones, local strength reduction is initiated respectively; by iteratively reducing the mechanical parameters of the starting zone and performing iterative calculations, the critical local reduction coefficient Fi corresponding to the overall instability of the system caused by each initial damage zone is obtained, and the minimum value of Fi is used as the basis for determining the main disaster source, and the main disaster source area is determined; (c) Taking the main disaster source as the core, the plastic zone expansion path caused by it is continuously tracked, and the newly formed plastic failure zone is included in the dynamic reduction range step by step until a continuous macroscopic failure surface is formed, thereby simulating a complete progressive failure process; (d) In steps (b) and (c), a partitioned progressive dynamic strength reduction strategy is introduced. The attenuation of rock mass strength parameters is controlled by field variables, and a functional relationship between the strength reduction coefficient Fv and the analysis step time is established. This makes the cohesion C and internal friction angle φ of the rock mass dynamically decrease as the analysis step progresses, until the system reaches the critical instability state. (e) During the iterative calculation process of steps (b) and (c) above, the multi-parameter fusion instability criterion is called in real time as the loop termination condition. The single non-convergence criterion is abandoned, and a multi-parameter fusion instability criterion based on the sudden change of characteristic point displacement, the plastic zone continuity rate and the sudden increase of system kinetic energy is established. The ultimate safety factor of the slope is determined by monitoring the inflection point of the data curve.

2. The stability analysis method for a goaf-slope system based on the dynamic strength reduction method according to claim 1, characterized in that, In step (a), the reduction factor used in the initial overall strength reduction calculation is slightly greater than 1.0, with a value of 1.05-1.2; the criteria for determining the initial damage area are: Under the aforementioned reduction factor, the equivalent plastic strain of the rock mass element reaches a preset threshold, or enters a shear or tensile yield state.

3. The stability analysis method for a goaf-slope system based on the dynamic strength reduction method according to claim 1, characterized in that, In step (b), the local strength reduction process is independent and exclusive; for each candidate initial damage zone, after resetting the model state, only the rock mechanics parameters within its region are reduced, including the tangent of cohesion c and internal friction angle φ. The specific steps for calculating the partial reduction are as follows: b1) Set the initial local reduction coefficient Fi = 0.8, the reduction step size ΔF, and the iteration step number i = 0; b2) calculating a reduction factor for the current iteration step i: F i = 0.8 + i x AF; b3) calculating the updated mechanical parameters of the rock mass in the initial damage zone according to the current reduction factor F i , the updated cohesion c i = c / F i , the updated internal friction angle φ i = arctan(tanφ / F i ), wherein c and φ are initial parameters of the rock mass; b4) Assign the updated parameters to the region cell and run the calculation until equilibrium is reached; b5) Determine whether the system is unstable based on the multi-parameter fusion instability criterion in step (e). If it is not unstable, let i = i + 1 and return to step b2) to continue the loop. If the system is unstable, the loop is terminated and the current F i as the critical local reduction factor for this region; The initial damage zone with the smallest critical local reduction coefficient is identified as the main source of disaster. When multiple initial damage zones have the same critical local reduction coefficient, the main disaster source is further determined by the expansion rate or expansion range of the plastic zone as an auxiliary criterion.

4. The stability analysis method for a goaf-slope system based on the dynamic strength reduction method according to claim 1, characterized in that, In step (c), the expansion of the dynamic reduction range follows the trigger-inclusion principle; In one calculation step, intensity reduction is performed only on the determined reduction range; After calculating the equilibrium, the cells that are adjacent to the current reduction zone and have newly entered the plastic state are identified to form a new plastic zone. In the next calculation step, this new plastic zone is merged into the dynamic reduction range, and the next level of strength reduction is applied to the updated entire reduction range. This process is repeated.

5. The goaf-slope system stability analysis method based on dynamic strength reduction method according to claim 1, characterized in that, In step (d), the formula for dynamic intensity reduction is as follows: where C i , φ i are the reduced cohesion and internal friction angle, C, φ are the initial parameters of rock mass, and Fv is the strength reduction factor; the dynamic increment is analyzed.

6. The goaf-slope system stability analysis method based on dynamic strength reduction method according to claim 1, characterized in that, In step (e), the criteria for determining the formation of the macroscopic failure surface are: the numerical calculation results do not converge, or the plastic zone extends from the toe of the slope or the boundary of the goaf to the top of the slope or other free surfaces, forming a continuous failure zone; the final stability coefficient of the system is determined by the critical local reduction coefficient of the main disaster source and the reduction in the subsequent progressive failure process.

7. The method of claim 6, wherein, In the multi-parameter fusion instability criterion, the weights of each criterion are fused through normalization to improve the stability and reliability of the instability determination.

8. Use of the method according to any one of claims 1 to 7, characterized in that, The method is applicable to the analysis of goaf-slope coupling systems composed of various lithologies and complex structural planes; it is also applicable to the evaluation of the long-term stability of the underlying goaf or the final slope of an open-pit mine affected by adjacent goaf activities.

Images