A method for predicting the risk of ice shelf disintegration and a monitoring system thereof

Abstract

This application relates to the field of polar marine disaster monitoring and early warning technology, and discloses a method and monitoring system for predicting ice shelf collapse risks. The method includes: establishing a sea-ice collaborative observation environment; obtaining the vertical displacement of the ice shelf using a base station-free carrier phase epoch difference algorithm; calculating wave propagation delay and establishing synchronous excitation response data pairs; dynamically correcting the equivalent bending stiffness of the ice shelf based on the energy transfer ratio evolution characteristics; combining real-time monitoring of the cross-sectional geometric curvature to invert and correct the instantaneous bending stress; statistically calculating the cumulative fatigue damage degree through stress cycles; constructing a ternary risk index system including instantaneous fracture, fatigue damage, and structural degradation; generating a comprehensive risk score and implementing graded control. This invention solves the problem of base station-free high-frequency deformation monitoring, corrects the stress calculation error caused by the constant stiffness assumption, and achieves accurate early warning of the entire process from damage accumulation to structural instability of the ice shelf.

Description

Technical Field

[0001] This invention relates to the field of polar marine disaster monitoring and early warning technology, specifically to a method for predicting ice shelf collapse risk and its monitoring system. Background Technology

[0002] Global warming has led to a significant increase in the instability of polar ice shelves. The frequent occurrence of ice shelf collapses not only affects sea-level changes but also poses a threat to polar scientific research and navigation safety. To assess the stability of ice shelves, it is usually necessary to monitor their dynamic response under wave loads in real time and to invert the internal stress state based on the monitoring data to predict the risk of collapse.

[0003] Global Navigation Satellite System (GNSS) technology is widely used in acquiring high-frequency deformation data of ice shelf surfaces. Traditional high-precision GNSS positioning (such as real-time dynamic differential (RTK)) heavily relies on reference stations deployed on stable bedrock. However, polar ice shelves are typically large-scale floating ice structures, with their edges often tens or even hundreds of kilometers away from inland bedrock, making it difficult to find absolutely stationary bedrock near the monitoring point to deploy reference stations. If the reference station is too far from the monitoring point, weakened spatial correlation leads to atmospheric errors that cannot be effectively eliminated, thus reducing positioning accuracy. This reliance on fixed base stations limits the application of existing monitoring technologies in the leading edge areas of floating ice shelves far from the coast, making it difficult to acquire high-precision dynamic bending deformation data.

[0004] In stress inversion and safety assessment based on monitoring data, existing technologies typically calculate internal bending stress by relying on the geometric curvature of the ice shelf surface. During the calculation process, most existing methods assume that the material properties of the ice shelf are constant, directly using fixed elastic modulus or bending stiffness values. However, under long-term wave cyclic loading, microcracks develop inside the ice shelf and gradually expand, leading to fatigue damage in the material, which macroscopically manifests as a gradual decrease in bending stiffness. If this degradation of material properties over time is ignored, and the constant stiffness parameters from the undamaged state are still used for inversion, the calculated stress values ​​will not accurately reflect the current load-bearing capacity of the ice shelf, resulting in distorted stress assessment results.

[0005] Existing ice shelf collapse early warning systems often focus on single-dimensional threshold judgments, such as monitoring only whether instantaneous stress or strain exceeds the critical fracture strength. While this assessment method can identify instantaneous fractures induced by extreme loads, it often overlooks the risks posed by the degradation of the structure's inherent performance. In reality, even if the instantaneous stress generated by external wave loads does not reach the theoretical fracture limit, the ice shelf may still collapse due to severe loss of structural stiffness caused by long-term cumulative damage. Current technologies lack a mechanism to identify this low-stress collapse mode induced by structural degradation, making it difficult to accurately predict the stability of the ice shelf throughout its entire life cycle. Summary of the Invention

[0006] To address the shortcomings of existing technologies, this invention provides a method for predicting ice shelf collapse risks and a monitoring system thereof. It solves the problems that existing polar ice shelf monitoring technologies are unable to achieve high-frequency deformation monitoring in the absence of stable bedrock reference stations, and that traditional stress inversion methods lead to assessment distortions because they ignore the characteristic that material stiffness decreases with damage accumulation. Furthermore, the single-dimensional risk criterion cannot identify low-stress collapse risks induced by structural degradation.

[0007] To achieve the above objectives, the present invention provides the following technical solution:

[0008] The first aspect of this invention provides a method for predicting the risk of ice shelf collapse, the method comprising the following steps:

[0009] Construct a multi-source collaborative observation environment. Deploy a BeiDou high-frequency positioning node array on the monitoring sections of the identified high-risk areas for ice shelf collapse, and deploy a wave monitoring array along the main propagation direction of the incident wave in the sea area outside the ice shelf leading edge to establish a synchronous observation environment for ocean waves and ice shelves.

[0010] Dynamic data of ocean waves and ice shelves were collected and processed. Wave surface displacement sequences of ocean waves were measured using a wave monitoring array; carrier phase observations were collected using a BeiDou high-frequency positioning node array; and the vertical displacement component sequence of the ice shelf surface was obtained by calculating the receiver position change vector between adjacent epochs and performing time integration without relying on a fixed land reference station, using a carrier phase epoch difference algorithm.

[0011] Spatiotemporal synchronization of wave data and ice shelf data was performed. Signal preprocessing was performed on the wave surface displacement sequence and the ice shelf vertical displacement component sequence. The dynamic time delay of wave propagation from the monitoring array to the ice shelf leading edge was calculated through cross-correlation analysis. Based on this dynamic time delay, time compensation was performed on the wave surface displacement sequence to establish a time-synchronized data pair of wave excitation and ice shelf response.

[0012] A dynamic correction model for the equivalent bending stiffness of the ice shelf was constructed. Based on synchronous data pairs, the excitation energy amplitude and response energy amplitude were extracted separately, and the real-time transfer ratio, reflecting the system's transfer characteristics, was calculated. According to the principles of structural dynamics, a nonlinear mapping relationship between the equivalent bending stiffness and the real-time transfer ratio was established. By utilizing the evolution trend of the real-time transfer ratio relative to the baseline state, the equivalent bending stiffness of the ice shelf was dynamically corrected to quantify the stiffness attenuation caused by crack propagation or cumulative damage within the ice shelf.

[0013] Inverted and corrected instantaneous bending stress. Using relative displacement data from adjacent nodes in the BeiDou high-frequency positioning node array, the instantaneous geometric curvature of the monitoring section is calculated using a differential method. Combining the dynamically corrected equivalent bending stiffness with the ice shelf thickness parameters, the instantaneous bending stress inside the ice shelf is inverted based on the elastic thin-plate bending theory, thereby correcting the stress calculation error caused by the assumption of constant material stiffness.

[0014] Calculate the cumulative fatigue damage of the ice shelf. Perform effective cyclic statistics on the instantaneous bending stress sequence obtained by inversion to identify the amplitude and number of stress cycles. Combine the fatigue life characteristics of ice materials with the linear cumulative damage theory to calculate the cumulative fatigue damage of the ice shelf under long-term wave cyclic loading.

[0015] A three-dimensional risk indicator system is constructed and a tiered early warning system is implemented. An instantaneous fracture risk indicator reflecting instantaneous strength, a fatigue damage risk indicator reflecting long-term cumulative effects, and a structural degradation risk indicator reflecting structural stiffness decay are calculated. A comprehensive risk score is generated through weighted fusion, and the comprehensive risk score and structural degradation risk indicator are compared with preset thresholds. Tiered control measures are implemented based on the comparison results.

[0016] In conjunction with the first aspect, in the first possible implementation, the process of obtaining the vertical displacement component sequence of the ice shelf using the carrier phase epoch difference algorithm specifically includes: constructing a carrier phase epoch difference observation equation with the receiver position change vector and receiver clock drift as unknown parameters; linearizing the observation equation and solving it using the least squares method to obtain the position change vector in the Earth-fixed coordinate system; projecting and transforming the position change vector to the station-centered horizontal coordinate system, extracting the vertical velocity component and performing cumulative integration to reconstruct the vertical displacement time series of the ice shelf surface.

[0017] In conjunction with the first aspect, in the second possible implementation, the process of constructing the dynamic correction model of the equivalent bending stiffness of the ice shelf specifically includes: setting an energy integration sliding window, calculating the root mean square value of the wave surface displacement sequence and the vertical displacement component sequence of the ice shelf within the window after time delay compensation, and defining their ratio as the real-time transfer ratio; establishing stiffness correction logic, setting the equivalent bending stiffness to be inversely proportional to the real-time transfer ratio, and when the real-time transfer ratio increases, correspondingly reducing the value of the equivalent bending stiffness, and ensuring that the corrected stiffness value does not exceed the initial theoretical bending stiffness.

[0018] In conjunction with the first aspect, in the third possible implementation method, the process of implementing hierarchical control specifically includes: presetting yellow warning thresholds, red warning thresholds, and stiffness critical thresholds; when the comprehensive risk score is between the yellow warning threshold and the red warning threshold and the structural degradation risk index does not exceed the limit, controlling the front-end monitoring equipment to increase the sampling frequency; when the comprehensive risk score exceeds the red warning threshold, or the structural degradation risk index exceeds the stiffness critical threshold, triggering a high-risk alarm, and predicting the drift trajectory of the collapsing iceberg based on the wave azimuth and ice shelf velocity.

[0019] A second aspect of the present invention provides an ice shelf collapse risk prediction and monitoring system, the monitoring system comprising a front-end sensing subsystem, a data transmission subsystem, and a data processing subsystem.

[0020] The front-end sensing subsystem includes a BeiDou high-frequency positioning node array deployed on the ice shelf monitoring section and a wave monitoring array deployed in the sea area at the front edge of the ice shelf, used to collect raw observation data.

[0021] The data transmission subsystem is used to establish a communication link between the front-end sensing subsystem and the data processing subsystem, and to transmit raw observation data and control commands.

[0022] The data processing subsystem is equipped with a computer program that, when executed by a processor, implements the following functional modules: a data calculation module, used to process carrier phase observations based on a carrier phase epoch difference algorithm to calculate the vertical displacement component sequence of the ice shelf, and to process wave observation data to calculate the wave surface displacement sequence; a spatiotemporal synchronization module, used to calculate the wave propagation delay and, accordingly, time-align the data to establish synchronized data pairs; a stiffness correction module, used to extract energy characteristics based on synchronized data pairs, calculate the real-time transfer ratio, and output the dynamically corrected equivalent bending stiffness; a stress inversion module, used to calculate the instantaneous geometric curvature of the monitoring section and, combined with the equivalent bending stiffness, invert the instantaneous bending stress; a damage assessment module, used to perform cyclic statistics on the stress sequence and calculate the cumulative fatigue damage degree; and a risk warning module, used to calculate a ternary risk index, generate a comprehensive risk score, and trigger corresponding graded control instructions according to the risk level.

[0023] This invention provides a method for predicting ice shelf collapse risk and a monitoring system thereof. It has the following beneficial effects:

[0024] 1. This invention utilizes a carrier phase epoch differential algorithm to achieve direct observation of high-frequency dynamic displacement on the ice shelf surface without relying on a fixed land-based reference station. This is achieved by solving the receiver position change vector between adjacent epochs and performing time integration. This technical feature effectively solves the engineering problem of the difficulty in setting up traditional differential positioning reference stations in polar floating ice shelf environments due to the lack of stable bedrock. It also eliminates the error accumulation caused by long baseline signal transmission and ensures that the monitoring system can accurately capture the minute bending deformation caused by wave excitation in a dynamic reference frame.

[0025] 2. This invention constructs a dynamic correction model for the equivalent bending stiffness of ice shelves. By calculating the transfer ratio between wave excitation energy and ice shelf response energy in real time, the equivalent bending stiffness parameters of the ice shelf are dynamically inverted and corrected. This method overcomes the limitation of traditional monitoring techniques that usually assume the elastic modulus of ice materials to be constant. It can quantify the material stiffness reduction caused by the propagation or cumulative damage of microcracks inside the ice shelf, thereby using stiffness values ​​that conform to the current actual state in the stress inversion process. This avoids stress calculation errors caused by overestimating the load-bearing capacity of the damaged ice shelf and significantly improves the accuracy of stress assessment.

[0026] 3. This invention establishes a three-dimensional risk index system that includes instantaneous fracture, fatigue damage, and structural degradation, and designs a stiffness-priority triggering logic with structural degradation as an independent criterion. This multi-dimensional assessment mechanism can not only provide early warning of the risk of instantaneous fracture caused by extreme wave loads, but also identify the potential collapse hazards induced by the severe loss of structural stiffness under low stress levels. It realizes comprehensive monitoring and early warning of the entire process of ice shelves from material damage accumulation to structural instability and collapse. Attached Figure Description

[0027] Figure 1 This is a flowchart of the method of the present invention;

[0028] Figure 2 This is a system architecture diagram of the present invention;

[0029] Figure 3 This is a schematic diagram of the wave excitation and stress response time series of the present invention;

[0030] Figure 4 This is a schematic diagram of the displacement and curvature calculation results of the present invention;

[0031] Figure 5 This is a schematic diagram illustrating the cumulative effect and trend analysis of the present invention;

[0032] Figure 6 This is a schematic diagram illustrating the risk scoring and level determination results of the present invention;

[0033] Figure 7 This is a schematic diagram of the normalized index and sub-score of the present invention. Detailed Implementation

[0034] The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0035] Please see the appendix Figure 1 - Appendix Figure 7 This invention provides a method for predicting the risk of ice shelf collapse and a monitoring system thereof. The monitoring system consists of a front-end sensing subsystem, a data transmission subsystem, and a data processing subsystem in its logical architecture.

[0036] The front-end sensing subsystem is configured to acquire external dynamic environment parameters and internal structural response parameters of the ice shelf. This subsystem includes a wave monitoring array and an ice surface deformation monitoring array. The wave monitoring array is deployed in open waters off the leading edge of the target ice shelf, upstream of the wave propagation path. Each wave monitoring array consists of at least two wave buoys spaced apart along the main propagation direction of the swell, used to measure wave height, period, and wave surface displacement sequences, providing external excitation input data for the data processing subsystem. The ice surface deformation monitoring array is deployed in key stress areas on the ice shelf surface and includes several BeiDou high-frequency positioning nodes. These nodes are linearly or arrayed along a direction perpendicular to the leading edge of the ice shelf, with the node spacing set based on the ice shelf thickness and theoretical bending wavelength to ensure the capture of the ice shelf's bending modes. Each BeiDou high-frequency positioning node integrates a Global Navigation Satellite System antenna, a BeiDou receiver, and a cryogenic power supply module for acquiring multi-frequency carrier phase observations.

[0037] The data transmission subsystem is configured to establish a bidirectional communication link between the front-end sensing subsystem and the data processing subsystem. Each wave buoy and BeiDou high-frequency positioning node transmits observation data frames to the aggregation gateway via a satellite short message communication module or a radio ad hoc network module. The data transmission subsystem uses BeiDou timing signals to uniformly timestamp all observation data to eliminate clock drift between different sensors and ensure strict synchronization of multi-source data on a time reference.

[0038] The data processing subsystem is configured to receive raw observation data, execute signal preprocessing, stiffness inversion, stress calculation, and risk assessment algorithms, and output early warning results. The data processing subsystem is a shore-based server or cloud computing platform with high concurrency processing capabilities, and it integrates data processing modules, spatiotemporal synchronization modules, stiffness correction modules, stress inversion modules, damage assessment modules, and risk early warning modules.

[0039] The ice shelf collapse risk prediction method described in this embodiment of the invention includes the following steps S100 to S700:

[0040] Step S100: Construct a multi-source collaborative observation environment. Based on historical remote sensing imagery and ice shelf thickness survey data, identify potential high-risk areas for ice shelf collapse at the ice shelf leading edge and determine the spatial locations of monitoring sections. Deploy ice surface deformation monitoring arrays on the selected monitoring sections, and simultaneously deploy wave monitoring arrays in the corresponding forward sea areas. During deployment, establish the spatial geometric relationship between the wave monitoring arrays and the ice shelf leading edge to facilitate subsequent calculations of the time delay of wave propagation to the ice shelf.

[0041] Step S200: Real-time acquisition and calculation of high-frequency data. The front-end sensing subsystem is activated, executing wave parameter acquisition and ice shelf deformation calculation in parallel. Wave buoys continuously acquire wave surface motion data at a preset sampling frequency, outputting sea state parameters including significant wave height, peak period, and wave surface displacement sequence. Simultaneously, BeiDou high-frequency positioning nodes acquire raw carrier phase observations. Utilizing a base station-free carrier phase epoch differential algorithm, without relying on fixed land-based base stations, the carrier phase observations of adjacent epochs are differentially processed to eliminate satellite clock errors and reduce tropospheric and ionospheric delay errors, enabling real-time calculation of the three-dimensional velocity vectors of each BeiDou high-frequency positioning node between adjacent epochs. Subsequently, the three-dimensional velocity vectors are integrated over time to obtain the raw three-dimensional displacement sequence of each BeiDou high-frequency positioning node in the Earth-fixed coordinate system, thereby achieving direct observation of the high-frequency dynamic motion of the floating ice shelf.

[0042] Step S300: Signal Preprocessing and Spatiotemporal Synchronization. The acquired raw data is cleaned and aligned. A high-pass filter is applied to the raw three-dimensional displacement sequence. By setting a cutoff frequency, the low-frequency background components caused by tidal changes, atmospheric pressure variations, and ice flow are separated, and the transient elastic response displacement components caused solely by wave excitation are extracted. Using BeiDou time as a unified reference, the wave propagation delay is calculated based on the physical distance between the wave monitoring array and the ice shelf leading edge, as well as the wave group velocity. The wave data is matched with the ice shelf response data corrected for wave propagation delay on the time axis to construct synchronized excitation and response data pairs, providing input for subsequent transfer function analysis.

[0043] Step S400: Dynamic Inversion of Ice Shelf Equivalent Stiffness Based on Load and Response Transfer Function. Using the synchronized wave excitation energy as input and the ice shelf deformation response energy as output, the real-time transfer ratio of the monitoring system is calculated. Based on thin plate bending theory, the deformation of the structure under a specific load is negatively correlated with the bending stiffness of the material. By monitoring the temporal evolution characteristics of the real-time transfer ratio, the current equivalent bending stiffness of the ice shelf is dynamically inverted. Specifically, a baseline state is set, which is either the initial stage of system initial operation or the stable period when the ice shelf shows no obvious signs of collapse. The transfer ratio and initial stiffness under the baseline state are obtained. When an abnormal increase in the response amplitude of the ice shelf is detected under the same wave excitation intensity, i.e., the real-time transfer ratio increases relative to the transfer ratio under the baseline state, it is determined that there is stiffness attenuation caused by crack propagation or cumulative damage inside the ice shelf. Based on this, the initial stiffness is corrected in real time, and the mechanical property model of the ice shelf is updated.

[0044] Step S500: Instantaneous Bending Stress Inversion Based on Dynamic Parameters. Combining geometric deformation observations with dynamic mechanical parameters, the true stress state inside the ice shelf is calculated. The instantaneous geometric curvature of the monitored section is calculated using the relative displacement difference between adjacent nodes in the ice surface deformation monitoring array. Subsequently, the corrected equivalent bending stiffness obtained in Step S400 is substituted into the bending constitutive equation to calculate the instantaneous bending stress inside the ice shelf. This step, by introducing a dynamically changing equivalent bending stiffness, corrects the stress calculation error caused by the assumption of a constant material modulus in traditional methods, thus improving the accuracy of stress state assessment for damaged ice shelves.

[0045] Step S600: Damage Accumulation and Fatigue Life Assessment. Based on the corrected instantaneous bending stress sequence, the long-term fatigue state of the ice shelf is assessed. A noise threshold is set to eliminate invalid micro-fluctuations, and effective stress cycles are extracted. The rainflow counting method is used to perform cyclic statistics on the instantaneous bending stress sequence to obtain different stress amplitudes and their corresponding cycle numbers. Combining the fatigue life curve of ice materials and the linear cumulative damage theory, the cumulative fatigue damage of the ice shelf under long-term wave alternating loads is calculated, quantifying the continuous impact of historical loads on structural stability.

[0046] Step S700: Multidimensional Collapse Risk Criterion and Graded Early Warning. A ternary joint criterion including instantaneous strength, fatigue damage, and stiffness degradation is constructed to comprehensively assess ice shelf stability. Instantaneous fracture risk index, fatigue damage risk index, and structural degradation risk index are calculated, and a comprehensive risk score is calculated based on the weight of each index's contribution to collapse. The comprehensive risk score and each individual index are compared with preset risk thresholds to determine the current risk level. Based on the determination result, a low-risk, medium-risk, or high-risk early warning signal is output, triggering corresponding system linkage operations. Specifically, when determined to be medium risk, the data processing subsystem issues an instruction to adjust the working mode of the BeiDou high-frequency positioning node to high-frequency sampling mode to obtain denser monitoring data; when determined to be high risk, the highest-level alarm signal is generated and the full data transmission mode is activated, while the potential iceberg drift trajectory is predicted based on wave and ice flow parameters.

[0047] Step S100 in this embodiment of the invention specifically involves the process of constructing a multi-source collaborative observation environment. The process of constructing a multi-source collaborative observation environment specifically includes sub-steps S110 to S130.

[0048] In sub-step S110, the monitoring section of the ice shelf and the key stress area are determined.

[0049] Multi-temporal synthetic aperture radar (SAR) and optical remote sensing images of the target ice shelf were acquired. Combined with ice shelf thickness survey data, the retreat rate of the ice shelf leading edge and the propagation rate of surface cracks were calculated. Areas where the ice shelf thickness reduction rate exceeded a preset threshold or where the crack distribution density was highest were selected as high-risk collapse zones. Within these high-risk zones, one or more monitoring sections perpendicular to the direction of the ice shelf leading edge were established. These monitoring sections covered the stress concentration zone from the ice shelf leading edge inland. The length of the monitoring sections was set according to the bending characteristic scale of the ice shelf to cover the main areas where elastic bending deformation of the ice shelf occurred.

[0050] In sub-step S120, an ice surface deformation monitoring array is deployed.

[0051] Several BeiDou high-frequency positioning nodes are installed along the monitoring section. These nodes are linearly distributed along the monitoring section, and the number of BeiDou high-frequency positioning nodes is... At least 3 ( This is to meet the mathematical requirements of calculating the second-order curvature using the three-point method in subsequent steps. The node spacing between two adjacent BeiDou high-frequency positioning nodes... Based on the estimated ice shelf bending wavelength Set the node spacing. The Nyquist sampling theorem must be satisfied, i.e. Preferred node spacing The range is from 100 meters to 500 meters. Each BeiDou high-frequency positioning node is fixed to the ice surface by an anchoring device and is equipped with a low-temperature adapted power module and a communication antenna.

[0052] In sub-step S130, a wave monitoring array is deployed. The wave monitoring array is deployed in open water outside the ice shelf leading edge. The array is positioned upstream of the wave propagation path to the monitoring section, and maintains a safe distance from the ice shelf leading edge. safe distance The preferred distance is 1 to 5 kilometers to prevent damage to the equipment in the wave monitoring array from breaking off icebergs. The wave monitoring array consists of at least two wave buoys, which are spaced apart along the main propagation direction of the incident wave according to their spatial relationship with the monitoring section, to measure the evolution characteristics of waves as they propagate to the ice shelf. The precise geographical coordinates of each wave buoy and each BeiDou high-frequency positioning node are determined using the Global Navigation Satellite System (GNSS), and the straight-line distance between the wave monitoring array and the ice surface deformation monitoring array is calculated as the value for calculating the wave propagation time delay in subsequent steps. The geometric input parameters.

[0053] In this embodiment of the invention, step S200 includes two parallel acquisition and calculation processes, wherein sub-step S210 specifically involves the real-time acquisition and processing of wave parameters. The wave parameter acquisition step further includes sub-steps S211 to S213.

[0054] In sub-step S211, raw data of ocean wave fluctuations are collected.

[0055] The wave buoy integrates an inertial measurement unit or accelerometer at a preset sampling rate. ( Continuous monitoring of the vertical acceleration of the sea surface The wave buoy's processing unit removes vertical acceleration. The gravitational acceleration component in the equation, and the vertical acceleration. A quadratic integration process is performed in the time domain. During this process, a high-pass filter is used to eliminate trend term drift caused by the accumulation of the integration constant, thereby obtaining time-series data on the vertical displacement of the sea surface over time. Simultaneously, the electronic compass or satellite compass built into the wave buoy synchronously acquires wave azimuth angle data.

[0056] In sub-step S212, sea state statistical parameters are calculated.

[0057] Fast Fourier Transform was performed on the time-series data of the collected sea surface vertical displacement over time to obtain the wave energy spectral density function. Based on wave energy spectral density function Calculate sea state characteristic statistics, including significant wave height. With peak period .

[0058] Significant wave height The calculation formula is:

[0059]

[0060] in, The zeroth moment of the wave spectrum, the zeroth moment The defining formula is:

[0061]

[0062] Peak period Defined as wave energy spectral density function The frequency corresponding to the maximum value The reciprocal of, that is Significant wave height With peak period The gating condition is used in the subsequent step S700 to determine whether the high-frequency sampling mode is triggered.

[0063] In sub-step S213, the synchronous wavefront displacement sequence is output.

[0064] The time-series data of sea surface vertical displacement over time obtained in sub-step S211 are resampled to align their time resolution with the system's unified time grid, thereby generating a wavefront displacement sequence. Using the BeiDou timing module mounted on a wave buoy, wave surface displacement sequences were generated. Each data point in the sequence is marked with a UTC timestamp. Wavefront displacement sequence with timestamps. As the external excitation input data of the system, it is sent to the data processing subsystem through the data transmission subsystem, and used in the subsequent step S400 for cross-correlation analysis and transfer function calculation with the ice shelf deformation response.

[0065] Step S220 in this embodiment of the invention specifically relates to the base station-free BeiDou high-frequency deformation calculation process. This base station-free BeiDou high-frequency deformation calculation process does not rely on fixed land-based reference stations, but instead utilizes the high correlation of a single BeiDou receiver in the time domain to obtain the high-frequency dynamic displacement of the ice shelf surface through carrier phase epoch differential technology. Step S220 further includes sub-steps S221 to S224.

[0066] In sub-step S221, each BeiDou high-frequency positioning node independently collects multi-frequency carrier phase observation values.

[0067] Set the sampling time interval to In the current epoch and the previous calendar year The Beidou receivers locked onto each other. visible satellites ( ). Obtain the BeiDou receiver's connection to the first Carrier phase observations of a satellite and The BeiDou receiver performs cycle slip detection and repair operations to ensure that integer ambiguity remains unchanged between epochs.

[0068] In sub-step S222, the carrier phase epoch differential observation equation is constructed and the receiver clock error and position change are solved.

[0069] Based on the inter-epoch difference of carrier phase observations within an extremely short time interval, integer ambiguity parameters are eliminated, and the residual effects of ionospheric and tropospheric delays are weakened. For the 1st... The carrier phase epoch differential observation equations are established for each satellite as follows:

[0070]

[0071] in, The carrier wavelength; For the calendar Carrier phase observations; For the receiver and the first The geometric distance between the satellites; The speed of light; For receiver clock bias; For satellite clock bias; This includes observation noise and unmodeled residual errors.

[0072] Calculate satellite positions and satellite clock biases using broadcast ephemeris or precise ephemeris To solve for the receiver position change, the carrier phase epoch differential observation equation is linearized and expanded. Satellite definition... The unit observation vector (direction cosine vector) to the receiver is Then the change in geometric distance can be approximated as:

[0073]

[0074] in, Let be the position change vector of the receiver in the Earth-fixed coordinate system; This represents the satellite's speed along the line of sight.

[0075] United The linearized equations of the visible satellites are constructed using the receiver position change vector. and receiver clock drift The system of error equations has unknown parameters. The least squares method is used to solve this system of error equations to obtain the position change vector of the BeiDou high-frequency positioning node in the Earth-fixed coordinate system. .

[0076] In sub-step S223, coordinate system transformation and velocity vector calculation are performed.

[0077] Obtain the approximate coordinates of the BeiDou high-frequency positioning node and calculate the transformation matrix from the Earth-fixed coordinate system to the station-centered horizontal coordinate system (ENU). The position change vector in the Earth-fixed coordinate system Projecting onto the station's horizontal coordinate system yields a displacement increment vector containing three components: east, north, and celestial. :

[0078]

[0079] Then calculate the instantaneous three-dimensional velocity vector of the current epoch. :

[0080]

[0081] in, The three components correspond to the eastward velocity, the northward velocity, and the vertical velocity, respectively.

[0082] In sub-step S224, the original three-dimensional displacement sequence is reconstructed by velocity integration.

[0083] To monitor the start time Based on this, the instantaneous three-dimensional velocity vectors of each epoch are accumulated and summed to obtain... The original relative displacement sequence at time 1 For vertical displacement components The integral formula for the vertical displacement component is:

[0084]

[0085] in, This represents the total number of epochs elapsed from the start time to the current time. , for Vertical velocity at any given moment.

[0086] The calculated original relative displacement sequence It contains high-frequency vibration information of the ice shelf and low-frequency trend terms caused by ice flow or tides, the original relative displacement sequence. This will serve as the input data for signal preprocessing in step S300. Through step S220, the ice shelf collapse risk prediction system achieves centimeter-level accuracy in monitoring relative deformation in floating ice shelf environments where stable base stations cannot be deployed.

[0087] Step S300 in this embodiment of the invention specifically relates to the construction and execution process of a wave-ice shelf spatiotemporal synchronization mechanism. The construction and execution process of the wave-ice shelf spatiotemporal synchronization mechanism specifically eliminates the differences in time reference and spatial propagation of multi-source observation data, constructing a synchronized dataset with physical causal correlation. Step S300 further includes sub-steps S310 to S320.

[0088] In sub-step S310, signal separation and high-pass filtering are performed on the original displacement sequence of the ice shelf.

[0089] The original relative displacement sequence (specifically the vertical displacement components) output in step S224 The data superimposed includes multiple frequency components, including high-frequency components of transient elastic bending caused by wave excitation and low-frequency trend components caused by ocean tides, atmospheric pressure loading effects, and the overall rheology of the ice shelf. The data processing subsystem uses the moving average method to calculate the low-frequency trend components. .

[0090] Define the sliding window length as ,time low-frequency trend components The calculation formula is:

[0091]

[0092] in, It is the integral variable.

[0093] Using vertical displacement components Subtract low-frequency trend components The vertical displacement components of the ice shelf after high-pass filtering were obtained. :

[0094]

[0095] Vertical displacement components of the ice shelf This only includes the purely elastic bending response of the ice shelf under leading-edge wave loads. The cutoff period of the high-pass filter. Set to 30 to 60 seconds.

[0096] In sub-step S320, spatiotemporal alignment and synchronization matching of wave data and ice shelf data are performed.

[0097] The spatiotemporal alignment and synchronization matching steps further include three processing stages: unified time reference, resampling alignment, and dynamic delay compensation.

[0098] Unified Time Reference and Resampling: Coordinated Universal Time (UTC) provided by the BeiDou Navigation Satellite System is used as the unified time reference for the entire system. A unified time raster sequence is defined. ,in The start time, For raster index, The standard sampling interval is (e.g., 1 second). Linear interpolation or cubic spline interpolation methods are used to interpolate the non-uniformly sampled wavefront displacement sequence. Vertical displacement component of the ice shelf Mapped to the same time raster The above yields the resampled wavefront displacement sequence. With ice shelf displacement sequence .

[0099] Wave propagation dynamic time delay calculation stage:

[0100] Due to the physical distance between the wave monitoring array and the leading edge of the ice shelf Furthermore, the wave propagation speed varies with sea state, and there is a time-varying propagation delay between the wave surface displacement sequence and the ice shelf displacement sequence. The data processing subsystem uses cross-correlation analysis for dynamic estimation. The maximum search latency is set to... Calculate the wavefront displacement sequence With ice shelf displacement sequence Different estimated delays Cross-correlation coefficients under :

[0101]

[0102] in, This refers to the number of data points within the correlation analysis window; and These are the mean values ​​of wave displacement and ice shelf displacement within the window, respectively. The test delay to be estimated has a range of values. .

[0103] The search results in cross-relationships. The time delay that reaches its maximum value is determined as the optimal propagation delay at the current moment. :

[0104]

[0105] Synchronous dataset construction phase:

[0106] Optimal propagation delay obtained through calculation Time-shift compensation is performed on the wavefront displacement sequence. Ice shelf displacement sequence at time and wavefront displacement sequence at time 1 Pairing is performed to form synchronous stimulus-response data pairs. The synchronous excitation-response data pair serves as the basis for the transfer function calculation and stiffness inversion in the subsequent step S400.

[0107] Step S400 in this embodiment of the invention further includes sub-step S410, which specifically involves the energy feature extraction process of excitation and response. The energy feature extraction process of excitation and response adopts the sliding time window integration method to extract the root mean square amplitude from the wavefront displacement sequence and the ice shelf vertical displacement component sequence after spatiotemporal synchronization.

[0108] In sub-step S410, the data processing subsystem constructs a sliding calculation window. The length of the energy integration sliding window is set to... Energy integration sliding window length The value should cover at least 3 to 5 complete significant wave cycles, and the preferred energy integral sliding window length is... The range is set to 60 to 300 seconds to eliminate the interference of random fluctuations in a single wave on the statistical results and ensure that the extracted energy features have statistical stability.

[0109] Based on the optimal propagation delay output in step S320 Compensated synchronous wavefront displacement sequence Calculate the current time excitation energy amplitude .

[0110] Excitation energy amplitude The calculation uses the root mean square algorithm, and the specific formula is as follows:

[0111]

[0112] in, For integration variables, For the integration time And after optimal propagation delay Corrected wavefront displacement values.

[0113] Based on the high-pass filtered vertical displacement components of the ice shelf output in step S310 Calculate the length of the sliding window for energy integration. Internal response energy amplitude .

[0114] Response energy amplitude The calculation formula is as follows:

[0115]

[0116] in, For the integration time The vertical displacement component of the ice shelf.

[0117] Calculated excitation energy amplitude With response energy amplitude This will be used as input data for calculating the real-time transfer ratio in the subsequent step S420.

[0118] In this embodiment of the invention, step S400 further includes sub-steps S420 and S430. The construction and solution process of the dynamic stiffness correction model is based on the evolution of the energy ratio between the input excitation and the output response, and the stiffness parameters in the subsequent stress inversion model are corrected.

[0119] In sub-step S420, the real-time transfer ratio is calculated.

[0120] The data processing subsystem utilizes the excitation energy amplitude extracted in step S410. With response energy amplitude Calculate the real-time transfer ratio Real-time transmission ratio Defined as the ratio of the response energy amplitude to the excitation energy amplitude. To prevent calculation divergence due to the excitation energy amplitude approaching zero, a regularization parameter is introduced into the denominator during the calculation. .

[0121] Real-time transmission ratio The calculation formula is:

[0122]

[0123] in, The preferred order of magnitude is a positive real number set according to the system noise level. to .

[0124] In sub-step S430, a dynamic correction model for the equivalent bending stiffness of the ice shelf is constructed and solved.

[0125] First, determine the baseline state parameters. Select a stable period during the initial startup of the monitoring system or a historically confirmed period of ice shelf health as the baseline time window. Calculate the baseline time window The average value of the internal real-time transfer ratio is used as the reference transfer ratio. Simultaneously, the initial theoretical bending stiffness was calculated based on the physical parameters of the ice shelf. .

[0126] Initial theoretical bending stiffness The calculation formula is:

[0127]

[0128] in, This is the elastic modulus of ice, typically ranging from 1 GPa to 9 GPa. The average thickness of the ice shelf at the monitoring section is obtained through ice radar surveys or historical data. The Poisson's ratio for ice is typically taken as 0.33.

[0129] Secondly, a nonlinear mapping relationship for stiffness attenuation is established.

[0130] Based on the principles of structural dynamics, the following dynamic correction equations are constructed:

[0131]

[0132] in, For a moment The equivalent bending stiffness; This is the stiffness sensitivity index.

[0133] Stiffness sensitivity index The stiffness sensitivity index is used to adjust the sensitivity of stiffness to changes in the transfer ratio. The range of values ​​is For ice shelf structures exhibiting significant brittleness, the preferred option is... .

[0134] Finally, the dynamic stiffness sequence is output.

[0135] The data processing subsystem calculates the equivalent bending stiffness. In the calculation process, if the equivalent bending stiffness The calculated value is greater than the initial theoretical bending stiffness. Then let Equivalent bending stiffness This will be used as a time-varying parameter for calculating the actual bending stress in step S500.

[0136] In this embodiment of the invention, step S500 further includes sub-steps S510 and S520, which use the ice surface geometric deformation observation data and dynamic stiffness parameters to calculate the corrected internal stress of the ice shelf based on the instantaneous bending stress inversion process of dynamic parameters.

[0137] In sub-step S510, the instantaneous geometric curvature of the monitored section is calculated.

[0138] The data processing subsystem selects BeiDou high-frequency positioning nodes at non-edge locations from the ice surface deformation monitoring array. As the calculation point, where , This represents the total number of nodes. (This refers to obtaining the number of BeiDou high-frequency positioning nodes.) and its preceding and following adjacent nodes and At the same time The vertical displacement component of the ice shelf (i.e., the output in step S310) ), respectively denoted as , and .

[0139] The three-point center difference method is used to calculate the BeiDou high-frequency positioning node. Instantaneous geometric curvature at point .

[0140] Instantaneous geometric curvature The calculation formula is:

[0141]

[0142] in, The physical distance between adjacent BeiDou high-frequency positioning nodes. This has been determined in step S120, the deployment phase.

[0143] In sub-step S520, the corrected instantaneous bending stress is inverted.

[0144] The time output in step S430 of the data processing subsystem call equivalent bending stiffness Combined with the instantaneous geometric curvature calculated in sub-step S510 The actual bending stress inside the ice shelf was calculated.

[0145] Based on the theory of elastic thin-plate bending, the maximum bending normal stress on the upper or lower surface of the ice shelf is calculated. The stress inversion formula is as follows:

[0146]

[0147] in, To monitor the average thickness of the ice shelf at the cross-section; This is the corrected instantaneous bending stress.

[0148] In the stress inversion formula, the equivalent bending stiffness These are dynamic parameters that vary over time. When cumulative damage to the ice shelf leads to an equivalent bending stiffness... As the numerical value decreases, for the same instantaneous geometric curvature The calculated instantaneous bending stress The numerical value decreases accordingly. The calculated instantaneous bending stress... This will be used as input data for fatigue life assessment in step S600.

[0149] Step S600 in this embodiment of the invention further includes sub-steps S610 and S620, specifically involving the cumulative calculation process of ice shelf fatigue damage. The cumulative calculation process of ice shelf fatigue damage utilizes the fatigue characteristics of ice materials to calculate the cumulative impact of long-term wave cyclic loads on the structural integrity of the ice shelf.

[0150] In sub-step S610, the effective stress cycles are extracted and the stress amplitude distribution is statistically analyzed.

[0151] The data processing subsystem receives the corrected instantaneous bending stress output in step S520. The data processing subsystem sets the stress threshold. stress threshold The value is taken as 1% to 5% of the critical fracture strength of the ice shelf. The data processing subsystem filters out instantaneous bending stress. The amplitude in the sequence is less than the stress threshold The fluctuation components are used to generate an effective stress sequence.

[0152] Subsequently, rainflow counting was used to perform cyclic statistics on the effective stress sequence. The data processing subsystem identified closed hysteresis loops in the effective stress sequence using rainflow counting and extracted stress cycles. Statistics were compiled from the monitoring start time. Up to the current moment The frequency of occurrence of stress amplitudes at various levels within a given time period. A stress amplitude spectrum is constructed to obtain... Different levels of stress amplitude and the corresponding cumulative number of loops ,in .

[0153] In sub-step S620, the cumulative fatigue damage is calculated.

[0154] The cumulative fatigue damage calculation process for ice shelves is based on Miner's linear cumulative damage theory and the SN fatigue curve of ice materials. First, the fatigue life equation of the ice materials is determined according to the physical properties of the target ice shelf (including ice temperature, salinity, and density).

[0155] The fatigue life equation is expressed as:

[0156]

[0157] in, For the first Level stress amplitude; To maintain a constant stress amplitude The limit number of cycles required for the ice shelf to undergo fatigue failure under certain conditions; This is the fatigue strength coefficient; It is the fatigue ductility index. and These are constants that were determined in advance through experiments on ice mechanics materials.

[0158] Next, calculate the current time. Cumulative fatigue damage Cumulative fatigue damage Defined as the sum of the damage ratios generated by each stress cycle, the calculation formula is as follows:

[0159]

[0160] in, For the statistical period Internal stress amplitude The actual cumulative number of cycles.

[0161] Cumulative fatigue damage It is a dimensionless scalar. When the cumulative fatigue damage... When the value is close to 1, the ice shelf is considered to be in a critical fatigue state. The calculated cumulative fatigue damage... As input parameters for the multidimensional risk criterion in step S700.

[0162] In this embodiment of the invention, step S700 further includes sub-step S710, which specifically relates to the construction process of a ternary risk indicator system. The construction process of the ternary risk indicator system includes sub-steps S711 to S713.

[0163] In sub-step S711, the instantaneous fracture risk index is calculated.

[0164] The data processing subsystem receives the corrected instantaneous bending stress output in step S520. The data processing subsystem obtains the critical fracture strength of the ice shelf. Critical fracture strength The ultimate stress value at which the ice material undergoes macroscopic brittle fracture, and the critical fracture strength. It is determined by in-situ ice mechanics test data or ice core laboratory test data from the relevant area.

[0165] Define instantaneous fracture risk index The ratio of the instantaneous bending stress amplitude to the critical fracture strength at the current moment is calculated using the following formula:

[0166]

[0167] in, This indicates the operation of taking the absolute value. Instantaneous fracture risk index. It is a dimensionless value. When the instantaneous fracture risk index... When the value exceeds 1.0, it is determined that the stress borne by the ice shelf at the current moment exceeds the theoretical strength of the material.

[0168] In sub-step S712, fatigue damage risk indicators are calculated.

[0169] The cumulative fatigue damage degree output in step S620 of the data processing subsystem call. Define fatigue injury risk indicators. Equal to cumulative fatigue damage The calculation formula is as follows:

[0170]

[0171] Fatigue injury risk indicators The numerical range is Fatigue injury risk indicators Used to quantify the cumulative damage state of ice shelves caused by the propagation of internal microcracks under historical wave cyclic loading.

[0172] In sub-step S713, the structural degradation risk index is calculated.

[0173] The data processing subsystem calls the equivalent bending stiffness calculated in step S430. Compared with the initial theoretical bending stiffness Define structural degradation risk indicators. The stiffness attenuation rate relative to the initial reference state is calculated using the following formula:

[0174]

[0175] According to the constraints in step S430 Structural degradation risk indicators The range of values ​​is When structural degradation risk indicators At that time, it was determined that the ice shelf stiffness had not degraded; when the structural degradation risk index... When the value approaches 1, the ice shelf is considered to have lost its ability to resist bending deformation. (Structural degradation risk index) As an independent physical state parameter, it is used to identify the risk of disintegration induced by structural material degradation at low stress levels.

[0176] In this embodiment of the invention, step S700 further includes sub-steps S720 and S730, which integrate multi-dimensional risk indicators into a quantitative score by combining the risk score calculation and graded early warning execution process, and automatically trigger graded responses based on the score results.

[0177] In sub-step S720, the comprehensive risk score is calculated.

[0178] The data processing subsystem employs a weighted fusion algorithm to process the instantaneous fracture risk index output in sub-step S710. Fatigue injury risk indicators and structural degradation risk indicators Perform comprehensive calculations. Define a comprehensive risk score. The calculation formula is as follows:

[0179]

[0180] in, These are the weight coefficients for the three indicators mentioned above, and the weight coefficients satisfy the normalization condition. In this embodiment, the weighting coefficients are set to have the following relationship: (For example, setting) This aims to increase the contribution of structural degradation risk indicators to the overall risk score.

[0181] To prevent sudden changes in scores caused by single transient noise, the data processing subsystem performs comprehensive risk scoring. Execution time smoothing and calculation of smoothed comprehensive score :

[0182]

[0183] in, To ensure a smooth time window length, a length of 30 to 60 seconds is preferred; It is the integral variable.

[0184] In sub-step S730, risk classification determination and system linkage are performed.

[0185] The data processing subsystem presets two risk classification thresholds: a yellow warning threshold. With red alert threshold And the numerical relationship satisfies Meanwhile, the data processing subsystem presets a critical stiffness threshold. The data processing subsystem will calculate the smoothed comprehensive score in real time. and structural degradation risk indicators The data is compared with the thresholds mentioned above, and a grading operation is performed based on the comparison results.

[0186] Level 1: Low-risk status (green).

[0187] When the condition is met and At this point, the ice shelf is determined to be in a stable period. The data processing subsystem controls the wave monitoring array and the ice surface deformation monitoring array to operate in a low-power mode, setting the data upload interval to a long period (e.g., 60 minutes) to reduce equipment power consumption.

[0188] Level 2: Medium risk (yellow).

[0189] When the condition is met and At that time, the ice shelf was determined to be in a period of risk accumulation. The data processing subsystem generates the first-level control command and sends it to each Beidou high-frequency positioning node and wave buoy through the data transmission subsystem. The first-level control command controls the front-end equipment to enter the high-frequency sampling mode, increasing the sampling frequency to 2 to 5 times the default frequency, and shortening the data upload interval to the minute level (e.g., 1 to 5 minutes).

[0190] Level 3: High-risk status (red).

[0191] When the condition is met or independent conditions If any one of these conditions is met, the ice shelf is determined to be in a critical disintegration phase. Among these, the independent condition... This constitutes a stiffness-priority triggering logic, which means that a red warning is forcibly triggered when the stress index has not exceeded the limit but the structural degradation risk index has exceeded the stiffness critical threshold.

[0192] In the third-level state, the data processing subsystem performs the following operations:

[0193] Generate the highest-level alarm signal and push it to the shore-based monitoring terminal via satellite link;

[0194] Enable full data storage and real-time transparent transmission mode, and disable front-end data compression and downsampling processing;

[0195] Based on the wave azimuth information obtained in step S210 and the instantaneous three-dimensional velocity vector of the ice shelf obtained in step S223 The kinematic extrapolation method was used to predict the initial drift trajectory of a potential iceberg that was breaking apart, and the predicted trajectory was superimposed on the electronic nautical chart interface for display.

[0196] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A method for predicting the risk of ice shelf collapse, characterized in that, Includes the following steps: Step S100: Construct a multi-source collaborative observation environment, deploy a BeiDou high-frequency positioning node array on the monitoring section of the high-risk area of ​​ice shelf collapse, and deploy a wave monitoring array along the main propagation direction of the incident wave in the sea area outside the ice shelf front. Step S200: Collect and calculate dynamic data of ocean waves and ice shelves, obtain wave surface displacement sequence through the wave monitoring array, and obtain the vertical displacement component sequence of ice shelf surface through the Beidou high-frequency positioning node array using carrier phase epoch difference algorithm. Step S300: Perform spatiotemporal synchronization of wave data and ice shelf data, calculate the dynamic time delay of wave propagation to the ice shelf, and establish a time delay compensated synchronous excitation response data pair; Step S400: Construct a dynamic correction model for the equivalent bending stiffness of the ice shelf, extract energy characteristics based on the synchronous excitation response data, calculate the real-time transfer ratio, and dynamically correct the equivalent bending stiffness of the ice shelf according to the evolution trend of the real-time transfer ratio. Step S500: Invert the corrected instantaneous bending stress, calculate the instantaneous geometric curvature using the monitoring data of the Beidou high-frequency positioning node array, and invert the instantaneous bending stress inside the ice shelf in combination with the equivalent bending stiffness. Step S600: Calculate the cumulative fatigue damage of the ice shelf by performing cyclic statistics on the instantaneous bending stress and calculating the cumulative fatigue damage based on the fatigue characteristics of the ice material; Step S700: Construct a three-element risk index system and implement graded early warning. Calculate the instantaneous fracture risk index, fatigue damage risk index and structural degradation risk index, integrate them to generate a comprehensive risk score, and implement graded control based on the comprehensive risk score and the structural degradation risk index.

2. The method for predicting ice shelf collapse risk according to claim 1, characterized in that, In step S200, the vertical displacement component sequence of the ice shelf surface is obtained by using the carrier phase epoch difference algorithm through the BeiDou high-frequency positioning node array, specifically including: Each of the aforementioned BeiDou high-frequency positioning nodes independently collects multi-frequency carrier phase observation values; a carrier phase epoch differential observation equation is constructed, wherein the carrier phase epoch differential observation equation uses the receiver position change vector and the receiver clock drift as unknown parameters; The carrier phase epoch differential observation equation is linearized, and the position change vector in the Earth-fixed coordinate system is solved using the least squares method. The position change vector is projected onto the station-centered horizontal coordinate system to obtain the instantaneous three-dimensional velocity vector containing the vertical component. The vertical velocity component in the instantaneous three-dimensional velocity vector is integrated over time to reconstruct the sequence of vertical displacement components of the ice shelf.

3. The method for predicting ice shelf collapse risk according to claim 1, characterized in that, The calculation of the dynamic time delay of wave propagation to the ice shelf in step S300 specifically includes: The wavefront displacement sequence and the ice shelf vertical displacement component sequence are subjected to detrending and high-pass filtering, and then uniformly resampled to the same time grid. Set a maximum search delay range, calculate the cross-correlation coefficient between the wavefront displacement sequence and the ice shelf vertical displacement component sequence under different estimated delays; search for the time delay corresponding to the cross-correlation coefficient reaching its maximum value, and determine the time delay as the dynamic delay.

4. The method for predicting ice shelf collapse risk according to claim 1, characterized in that, The step S400, which involves constructing a dynamic correction model for the equivalent bending stiffness of the ice shelf, specifically includes: Set the length of the energy integration sliding window, and perform root mean square calculations on the wavefront displacement sequence and the ice shelf vertical displacement component sequence after time delay compensation within the length of the energy integration sliding window to obtain the excitation energy amplitude and the response energy amplitude. The ratio of the response energy amplitude to the excitation energy amplitude is calculated to obtain the real-time transfer ratio; Obtain the baseline transfer ratio and initial theoretical bending stiffness within the baseline time window; A nonlinear mapping relationship for stiffness attenuation is established, wherein the equivalent bending stiffness is inversely proportional to the real-time transfer ratio. When the real-time transfer ratio increases relative to the reference transfer ratio, the value of the equivalent bending stiffness is reduced.

5. The method for predicting ice shelf collapse risk according to claim 1, characterized in that, The inverted and corrected instantaneous bending stress in step S500 specifically includes: Select the target BeiDou high-frequency positioning node and its adjacent nodes on the monitoring section; use the three-point center difference method to calculate the instantaneous geometric curvature at the target BeiDou high-frequency positioning node based on the vertical displacement component values ​​of the ice shelf at the target BeiDou high-frequency positioning node and its adjacent nodes; based on the elastic thin plate bending theory, calculate the product of the equivalent bending stiffness and the instantaneous geometric curvature, and combine it with the ice shelf thickness parameter to deduce the corrected instantaneous bending stress.

6. The method for predicting ice shelf collapse risk according to claim 1, characterized in that, The calculation of cumulative fatigue damage to the ice shelf in step S600 specifically includes: A stress threshold is set, and the fluctuation components with amplitudes smaller than the stress threshold in the instantaneous bending stress are filtered out to generate an effective stress sequence. The effective stress sequence is statistically analyzed using the rainflow counting method to identify closed hysteresis loops and obtain the stress amplitudes at each level and their corresponding cumulative cycle counts. Based on the fatigue life equation of ice materials, the limit number of cycles corresponding to each stress amplitude is calculated; using the linear cumulative damage theory, the ratio of the actual cumulative number of cycles to the limit number of cycles for each stress amplitude is summed to obtain the cumulative fatigue damage degree.

7. The method for predicting ice shelf collapse risk according to claim 1, characterized in that, The construction of the ternary risk indicator system in step S700 specifically includes: Calculate the instantaneous fracture risk index, which is the ratio of the amplitude of the instantaneous bending stress to the critical fracture strength of the ice shelf; Calculate the fatigue damage risk index, where the value of the fatigue damage risk index is equal to the cumulative fatigue damage degree. Calculate the structural degradation risk index, which is the attenuation ratio of the equivalent bending stiffness to the initial theoretical bending stiffness.

8. The method for predicting ice shelf collapse risk according to claim 7, characterized in that, The step S700, which involves generating a comprehensive risk score, specifically includes: Weighting coefficients are assigned to the instantaneous fracture risk index, the fatigue damage risk index, and the structural degradation risk index, respectively; the weighted summation result is smoothed by moving average to generate the comprehensive risk score; wherein the weighting coefficient corresponding to the structural degradation risk index is greater than the weighting coefficients corresponding to the instantaneous fracture risk index and the fatigue damage risk index.

9. The method for predicting ice shelf collapse risk according to claim 8, characterized in that, The hierarchical control implemented in step S700 specifically includes: Preset yellow warning threshold, red warning threshold, and stiffness critical threshold; When the comprehensive risk score is between the yellow warning threshold and the red warning threshold and the structural degradation risk index does not exceed the stiffness critical threshold, a medium-risk state is triggered, and the Beidou high-frequency positioning node array is controlled to enter the high-frequency sampling mode. When the comprehensive risk score exceeds the red warning threshold, or the structural degradation risk index exceeds the stiffness critical threshold, a high-risk state is triggered, an alarm signal is generated, and the drift trajectory of the collapsing iceberg is predicted based on the wave azimuth and ice shelf velocity.

10. A system for predicting and monitoring the risk of ice shelf collapse, used to implement the method for predicting the risk of ice shelf collapse as described in any one of claims 1 to 9, characterized in that, The monitoring system includes: The front-end sensing subsystem includes a BeiDou high-frequency positioning node array deployed on the ice shelf monitoring section and a wave monitoring array deployed in the sea area at the edge of the ice shelf, which are used to collect raw observation data. A data transmission subsystem is used to transmit the raw observation data to the data processing subsystem; The data processing subsystem is equipped with a computer program. When executed by a processor, the computer program performs the functions of the following modules: a data calculation module, used to calculate the vertical displacement component sequence of the ice shelf based on the carrier phase epoch difference algorithm, and to calculate the wave surface displacement sequence; a spatiotemporal synchronization module, used to calculate the wave propagation delay and establish a synchronous excitation response data pair; a stiffness correction module, used to calculate the real-time transfer ratio based on energy feature extraction and output the equivalent bending stiffness; a stress inversion module, used to calculate the instantaneous geometric curvature and invert the instantaneous bending stress in combination with the equivalent bending stiffness; a damage assessment module, used to calculate the cumulative fatigue damage degree; and a risk warning module, used to calculate the ternary risk index, generate a comprehensive risk score, and trigger graded control instructions.

Images