Tide and wave joint test method and system based on multi-source data in port engineering construction

Abstract

The invention relates to the technical field of port engineering construction, in particular to a tide wave joint test method and system based on multi-source data in port engineering construction. At present, the problems of insufficient data utilization, low test precision and the like exist in tide and wave joint test in port engineering construction. The method comprises the following steps: constructing a nested unstructured grid and dynamically configuring discrete parameters of a wave direction spectrum; a variation assimilation algorithm is adopted to optimize hydrodynamic parameters; solving a wave flow coupling equation to generate a combined field data set; the pile foundation stress is calculated in combination with an improved Morison equation, and a safety evaluation report is output through basin resonance characteristic analysis; the system sets corresponding modules to realize various functions based on the method. The device is mainly used for precise combined testing of tidal waves in port engineering construction so as to guarantee the safety of an engineering structure.

Description

Technical Field

[0001] The present invention relates to the technical field of port engineering construction, and more particularly to a method and system for joint testing of tidal currents and waves based on multi-source data in port engineering construction. Background Art

[0002] During port construction, accurately understanding the characteristics of tidal currents and waves is crucial to the project's stability, safety, and long-term operation. Tidal currents and waves not only affect port navigation conditions but also exert complex dynamic effects on port hydraulic structures such as docks, breakwaters, and pile foundations.

[0003] Traditional tidal wave testing methods have many limitations. On the one hand, the means of data acquisition are relatively simple, and often only limited local data can be obtained, making it difficult to fully reflect the complex hydrodynamic environment of the entire port area. For example, relying on only a few observation points to collect flow velocity and wave height data cannot accurately capture the spatial variation characteristics of flow velocity and waves in the area, resulting in an inaccurate description of the overall hydrodynamic field. On the other hand, there are difficulties in the fusion processing of different types of data. Multi-source data such as seabed topography data, coastline data, and hydrological and meteorological data are difficult to effectively integrate and analyze due to their different sources and inconsistent temporal and spatial resolutions. This makes it impossible to fully utilize all aspects of information when constructing a comprehensive hydrodynamic model, affecting the accuracy and reliability of the model.

[0004] When constructing numerical models to simulate tidal waves, existing methods lack precision in meshing and parameter setting. Traditional meshing methods cannot adapt well to complex seabed topography and shoreline features. In areas with drastic topography changes or in critical areas, such as around pile foundations, they lack sufficient mesh resolution, resulting in insufficient computational accuracy. Furthermore, the discrete settings for key parameters, such as the wave directional spectrum, lack specificity and cannot be properly adjusted based on actual hydrodynamic conditions and energy distribution, leading to deviations between simulation results and actual conditions.

[0005] In addition, when considering the wave-current coupling effect, the existing calculation methods fail to fully consider the complex mutual influence mechanism between the two. For example, when calculating the radiation stress gradient correction term, it cannot accurately reflect the impact of factors such as wave breaking and terrain changes on the tidal field, resulting in inaccurate simulation of the wave-current combined field, which in turn affects the accurate assessment of the stress on the port engineering structure. Moreover, the analysis methods for special phenomena such as harbor resonance are not perfect, making it difficult to accurately predict the resonance risk and unable to provide strong support for the optimized design of port engineering. These problems seriously restrict the quality and safety of port engineering construction, and a new tidal wave combined testing method is urgently needed to solve the above problems. Summary of the Invention

[0006] An object of the present invention is to solve at least the above problems and to provide at least the advantages which will be described hereinafter.

[0007] In order to achieve these purposes and other advantages according to the present invention, a method for joint testing of tidal currents and waves in port engineering construction based on multi-source data is provided, comprising the following steps:

[0008] S1. Acquire seabed topography data and shoreline point cloud data, collect multi-layer wind speed and direction data, current velocity profile data, and wave direction spectrum data;

[0009] S2, spatially registering the terrain data with the shoreline point cloud data to generate a seabed digital elevation model;

[0010] S3, perform spatiotemporal alignment processing on the multi-layer wind speed and direction data, velocity profile data, and wave direction spectrum data collected in S1 to generate a quality verification data set with temporal and spatial consistency;

[0011] S4, constructing an unstructured grid based on the seabed digital elevation model generated in S2, generating a nested spatial discrete model, and setting wave directional spectrum discrete parameters in the nested spatial discrete model;

[0012] S5, using the variational assimilation algorithm to optimize the quality verification data set of S3 and output the optimized hydrodynamic parameter combination;

[0013] S6, based on the wave directional spectrum discrete parameters of S4 and the hydrodynamic parameters optimized by S5, solve the wave-current coupling equation and calculate the radiation stress gradient correction term to generate the wave-current joint field data set;

[0014] S7. Configure the wave incident boundary conditions according to the wave direction spectrum data, synthesize the tide level boundary conditions according to the tidal harmonic constant, and combine the tide level and the design wave height to generate the tide wave combination condition;

[0015] S8, calculate the pile foundation stress parameters based on the wave-current combined field data set generated in S6 and the tide-wave combined working condition in S7;

[0016] S9, based on the wave-current joint field data set generated by S6 and the tide-wave combined working condition of S7, uses the resonance characteristic equation to calculate the frequency matching degree and output the resonance amplification factor;

[0017] S10. Generate an engineering structure safety assessment report based on pile foundation stress parameters and resonance amplification factors.

[0018] Preferably, the step S4 specifically includes:

[0019] S4.1. Extract the terrain gradient field based on the seafloor digital elevation model, divide the grid density control area according to the gradient threshold, and set the core density area, transition area, and background sparse area;

[0020] S4.2. Generate an initial unstructured grid using the Delaunay triangulation algorithm. For the core densification area specified in S4.1, implement local secondary densification using the superimposed frontier advancing method, ensuring that the grid side lengths meet the following requirements: core area ≤ 2m, transition area ≤ 5m, and background area ≤ 10m.

[0021] S4.3. Construct a three-level nested hierarchy based on the unstructured grid generated in S4.2: a first-level global grid covering the entire computational domain with a resolution of 10-20m; a second-level local grid covering the port waters with a resolution of 5-8m; and a third-level refined grid focusing on the area surrounding the pile foundation cluster with a resolution of 1-3m.

[0022] S4.4, coupling the nesting boundaries of the three nesting levels in S4.3 by a bidirectional flux interpolation method, and using virtual units to transfer the gradient flux of hydrodynamic parameters, wherein the hydrodynamic parameters include flow velocity and wave height;

[0023] S4.5. Based on the hydrodynamic parameter combination optimized in S5, dynamic adaptive encryption is triggered in the velocity shear layer and wave breaking zone, and the grid topology is adjusted in real time according to the grid hierarchy structure in S4.3.

[0024] Preferably, the step S4 further includes:

[0025] S4.6. Extract the shoreline main axis direction θ based on the seabed digital elevation model from S2, using the range of θ ± 30° as the wave incident main direction interval. Based on the multi-layer wind speed and direction data, current velocity profile data, and wave direction spectrum data collected from S1, discretize the direction spectrum at intervals of 5°-10° within the wave incident main direction interval, and at intervals of 15°-20° in other areas.

[0026] S4.7. Based on the wave directional spectrum energy distribution collected in S1, three frequency bands are divided: the core frequency band frequency interval Δf1 = 0.02 Hz, the transition frequency band Δf2 = 0.05 Hz, and the edge frequency band Δf3 = 0.1 Hz;

[0027] S4.8. Combine the three-level nested hierarchy constructed in S4.3 and configure the discrete parameters hierarchically in the three-level nested grids: first-level global grid, directional spacing of 15°, frequency spacing Δf3; second-level local grid, directional spacing of 10°, frequency spacing Δf2; third-level refined grid, directional spacing of 5°, frequency spacing Δf1;

[0028] S4.9, introduces a wavenumber-water depth coupling constraint. Based on the seafloor digital elevation model in S2, when the local water depth h and wavelength L of the seafloor digital elevation model satisfy h / L≤0.05, trigger the diffraction enhancement mode of the directional spectrum, and improve the directional resolution of the area to 2°;

[0029] S4.10. Based on the velocity field assimilated in S5, the multi-directional wave interference compensation algorithm is activated for the area where the velocity-wave angle is greater than 60°, and the reverse secondary wave component is superimposed on the original directional spectrum discretized in S4.6.

[0030] Preferably, the step S5 specifically includes:

[0031] S5.1. Construct a two-stage variational assimilation framework: In the tidal assimilation stage, the bottom friction coefficient and vertical eddy viscosity coefficient are optimized using the velocity profile data in the quality verification dataset generated by S3 as constraints; in the wave assimilation stage, the JONSWAP spectral peak enhancement factor and breakup index are optimized using the wave directional spectrum data in the quality verification dataset generated by S3 as constraints;

[0032] S5.2. Set a terrain correction cost function, which is a weighted sum of three terms: the mean square error term between the observed flow velocity and the model-calculated flow velocity, a coupling constraint term for the seafloor topography gradient and bottom shear stress, and a parameter sparsity regularization term. The weight coefficient of the topography gradient term is 0.3 to 0.6 times the weight of the observation term.

[0033] S5.3. Use the quasi-Newton-L-BFGS algorithm to solve the cost function in S5.2. When the gradient descent rate for three consecutive iterations is less than the preset threshold, the step size is automatically reduced to 1 / 4 of the original value;

[0034] S5.4. For areas where velocity and wave data conflict, dynamic weights are applied according to the following rules: in the time dimension, the weight decays exponentially with a half-life of 2-3 hours, centered on the conflict start time; in the spatial dimension, the weight decays to 25%-30% of the center area within a 50-60 meter radius centered on the pile foundation.

[0035] S5.5. Output the optimized parameter combination and confidence interval. When the coefficient of variation of tidal parameters is greater than 15% or the coefficient of variation of wave parameters is greater than 25%, manual review is triggered.

[0036] Preferably, the step S6 specifically includes:

[0037] S6.1. Solve the wave-current coupling equations in stages: Based on the discrete parameters of the wave directional spectrum from S4, in the wave-dominated phase, the gentle slope equation is used to calculate the wave field, outputting wave height, wave direction, and breaking energy distribution. Based on the hydrodynamic parameters optimized in S5, in the tidal-dominated phase, the three-dimensional shallow water equations are solved to obtain the velocity and water level fields. In the joint correction phase, the wave radiation stress gradient is injected as a source term into the tidal current equation, and the velocity field is fed back into the wave refraction calculation.

[0038] S6.2. Calculate the radiation stress gradient correction term: Use a turbulent kinetic energy compensation model in the wave breaking zone to convert 10%-15% of the dissipated energy from breaking into radiation stress increments. Dynamically adjust the correction term calculation step size based on the terrain gradient of the seafloor digital elevation model in S2. The step size should be ≤2m in steep slopes and ≤10m in flat areas.

[0039] S6.3. Perform bidirectional coupling iterations: After each wave-current coupling calculation in S6.1, compare the relative error of wave height and velocity deviation between two consecutive iterations. Terminate the calculation if the threshold is not exceeded for three consecutive iterations. Otherwise, update the radiation stress field and iterate again based on the radiation stress gradient correction term in S6.2.

[0040] S6.4. Generate a wave-current joint field dataset, which includes: a spatial vector field, including flow velocity, flow direction, wave height, and wave direction; an energy scalar field, including radiation stress intensity and turbulent kinetic energy density; and a dynamic characteristic field, including the location of the breaking zone and the vortex core area.

[0041] Preferably, the pile foundation force calculation in step S8 specifically includes:

[0042] S8.1. Input the wave-current combined field data set generated in S6 into the improved Morison equation at a time step of 0.1 seconds. Dynamically adjust the resistance coefficient C according to the angle between the flow direction and the pile axis and the turbulent kinetic energy in the wave-current combined field data of S6. d and the inertia coefficient C m :

[0043] When the angle between the flow direction and the pile axis is greater than 45°, the resistance coefficient C is dynamically adjusted. d Increase by 10% to 15%; in the vortex shedding area, the inertia coefficient C m Reduce by 5% to 8%, and the vortex shedding zone is turbulent kinetic energy>0.3J / m 3 ;

[0044] S8.2. Using the coefficients adjusted in S8.1, calculate the shear force and bending moment time history curves for each section of the pile foundation, and then extract the maximum value, standard deviation, and fatigue damage index.

[0045] Preferably, the step S9 specifically includes:

[0046] S9.1. Establish a three-dimensional characteristic equation based on the seabed digital elevation model generated in S2 and the harbor contour line, and solve the first six resonance modes with a frequency range of 0.01 to 0.5 Hz;

[0047] S9.2. Compare the dominant frequency of the wave spectrum in the wave-current joint field dataset generated in S6 with the resonant modal frequency calculated in S9.1 and calculate the degree of match: when the dominant frequency difference is ≤0.02 Hz, it is determined to be a strong resonance risk, and the amplification factor is set to 1.5 to 3.0; when the dominant frequency difference is >0.1 Hz, it is determined to be no resonance, and the amplification factor is set to 1.0;

[0048] S9.3. Determine the location of the breakwater opening based on the seafloor digital elevation model generated in S2. Add a 10% to 20% increment to the amplification factor obtained in S9.2 to account for the effects of secondary reflections superimposed on this area.

[0049] A tidal wave combined testing system based on the above method is provided, comprising:

[0050] A data acquisition module is used to obtain seabed topography data, shoreline point cloud data, multi-layer wind speed and direction data, flow velocity profile data, and wave direction spectrum data;

[0051] The elevation model generation module is used to spatially register terrain data with shoreline point cloud data to generate a seabed digital elevation model;

[0052] A data processing module is used to perform spatiotemporal alignment processing on the collected multi-layer wind speed and direction data, flow velocity profile data, and wave direction spectrum data to generate a quality verification data set with temporal and spatial consistency;

[0053] A model construction module is used to construct an unstructured grid based on the seabed digital elevation model, generate a nested spatial discrete model, and set wave directional spectrum discrete parameters in the nested spatial discrete model;

[0054] Parameter optimization module, which is used to optimize the S3 quality verification data set using the variational assimilation algorithm and output the optimized hydrodynamic parameter combination;

[0055] The wave-current combined field generation module is used to solve the wave-current coupling equation and calculate the radiation stress gradient correction term based on the combination of wave directional spectrum discrete parameters and hydrodynamic parameters to generate a wave-current combined field data set;

[0056] The working condition construction module is used to configure the wave incident boundary conditions according to the wave direction spectrum data, synthesize the tide level boundary conditions according to the tidal harmonic constant, and combine the tide level with the design wave height to generate the tidal wave combination working condition;

[0057] The pile foundation force calculation module is used to calculate the pile foundation force parameters based on the generated wave-current joint field data set and the tidal wave combined working condition;

[0058] The resonance factor calculation module is used to calculate the frequency matching degree using the resonance characteristic equation based on the generated wave-current joint field data set and the tide-wave combination condition, and output the resonance amplification factor;

[0059] The assessment report generation module is used to generate an engineering structure safety assessment report based on pile foundation stress parameters and resonance amplification factors.

[0060] The present invention has at least the following beneficial effects:

[0061] First, the use of unstructured grid nesting modeling and the L-BFGS parameter assimilation algorithm increases the grid generation speed by 50% and shortens the convergence time of hydrodynamic parameter optimization to 40% of traditional methods, meeting the real-time simulation requirements under complex working conditions.

[0062] Second, based on the improved Morison equation and resonant frequency matching algorithm, the pile foundation force prediction error is ≤8%, and the resonance risk identification accuracy reaches 95%, providing quantitative safety thresholds for port structure design (such as a bending moment ratio threshold of 1.0 and a resonance factor threshold of 0.8), reducing the risk of structural failure.

[0063] Third, through radiation stress gradient correction and bidirectional iterative coupling algorithm, the contribution of turbulent kinetic energy in the wave breaking zone to tidal kinetic energy can be quantified (compensating 10%-15% of dissipated energy), realizing dynamic interactive simulation of velocity field and wave field. Compared with the traditional one-way coupling method, the comprehensive error of pile foundation force calculation is reduced by 12%-18%.

[0064] Fourth, the combined tide and wave conditions constructed based on the asynchronous extreme value method can simulate extreme scenarios where the phase difference between the rising moment and the incoming wave crest is greater than 15 minutes. The combined test system's comprehensive response accuracy to the breakwater diffraction effect (direction fine-tuning 2°-3°) and the bay top terrain attenuation (amplitude attenuation 10%-20%) is improved by 25%, supporting engineering verification of extreme conditions that occur once every 50 years.

[0065] Fifth, the multi-scale dynamic process is analyzed collaboratively, using a three-level nested grid and directional spectrum graded discretization to capture vortex shedding (turbulent kinetic energy > 0.3 J / m) at the local pile foundation (1-3 m resolution). 3 ) details while maintaining tidal phase synchronization in the global domain (10-20m resolution), achieving multi-scale joint simulation from millimeter-level wave breaking to kilometer-level tidal propagation, and reducing computing resource consumption by 40% compared to traditional uniform grids.

[0066] Other advantages, objectives and features of the present invention will be reflected in part through the following description, and in part will be understood by those skilled in the art through study and practice of the present invention. BRIEF DESCRIPTION OF THE DRAWINGS

[0067] Figure 1 This is a framework flow chart of the joint testing system according to one of the technical solutions of the present invention. DETAILED DESCRIPTION

[0068] The present invention will be described in further detail below in conjunction with the accompanying drawings so that those skilled in the art can implement the invention with reference to the description.

[0069] It should be noted that the experimental methods described in the following embodiments are conventional methods unless otherwise specified, and the reagents and materials are commercially available unless otherwise specified; in the description of the present invention, the orientation or positional relationship indicated by the terms is based on the orientation or positional relationship shown in the accompanying drawings, which is only for the convenience of describing the present invention and simplifying the description, and does not indicate or imply that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and therefore cannot be understood as a limitation on the present invention.

[0070] like Figure 1 As shown, the present invention provides a tidal wave joint testing method based on multi-source data in port engineering construction, comprising the following steps:

[0071] S1. Acquire seabed topography data and shoreline point cloud data, collect multi-layer wind speed and direction data, current velocity profile data, and wave directional spectrum data. Specifically, during the data acquisition phase, three-dimensional shoreline point cloud data is acquired using a lidar, and seabed topography data is acquired using a multi-beam bathymetry system. The following time series data are simultaneously collected: a) multi-layer gradient wind speed and direction data (10m / 30m / 50m altitude layers); b) ADCP current velocity profile data (vertical resolution 0.5m); and c) directional spectrum data measured by a wave buoy (frequency resolution 0.01Hz, directional resolution 5°).

[0072] S2. spatially registering the terrain data with the shoreline point cloud data to generate a seabed digital elevation model. Specifically, terrain fusion modeling: spatially registering the shoreline 3D point cloud data with the seabed terrain data, using the Kriging interpolation algorithm to generate a seabed digital elevation model with a resolution of 0.5 m × 0.5 m, and extracting the terrain gradient field.

[0073] S3. Perform spatiotemporal alignment on the multi-layer wind speed and direction data, velocity profile data, and wave direction spectrum data collected in S1 to generate a quality verification dataset with temporal and spatial consistency. Specifically, the spatiotemporal alignment of multi-source data involves: a) performing clock synchronization based on UTC timestamps, unifying the sampling interval to 10 minutes; b) using the Kalman filter algorithm to correct the phase lag of wind speed and wave data; and c) filling in missing data segments in the velocity profile using cubic spline interpolation to generate a quality verification dataset with temporal and spatial consistency.

[0074] S4. Construct an unstructured grid based on the seabed digital elevation model generated by S2, generate a nested spatial discrete model, and set wave directional spectrum discrete parameters in the nested spatial discrete model. Specifically, based on the seabed digital elevation model generated by S2, perform the following: a) divide the core densification area (gradient > 5%), transition area (1%-5%), and background area (<1%) according to the terrain gradient; b) use Delaunay triangulation to construct a three-level nested unstructured grid: global grid (resolution 20m), port operation area grid (resolution 8m), and pile foundation group densification grid (resolution 2m); c) according to the wave directional spectrum data characteristics of S1, set the spectrum discrete parameters: main wave direction interval directional resolution 5°, secondary wave direction interval directional resolution 15°, frequency resolution hierarchical configuration (core frequency band 0.02Hz, edge frequency band 0.1Hz);

[0075] S5. Use the variational assimilation algorithm to optimize the quality verification data set of S3 and output the optimized hydrodynamic parameter combination. Specifically, the data assimilation optimization: a) optimize the bottom friction coefficient and eddy viscosity coefficient using the ADCP velocity profile as a constraint; b) optimize the JONSWAP spectrum peak factor using the wave direction spectrum as a constraint; c) correct the wind stress parameter by combining the gradient wind speed data and output the optimized hydrodynamic parameter combination and confidence interval;

[0076] S6: Based on the wave directional spectrum discrete parameters from S4 and the hydrodynamic parameters optimized in S5, the wave-current coupling equation is solved and the radiation stress gradient correction term is calculated to generate a wave-current joint field dataset. Specifically, the wave-current coupling calculation: a) solves the shallow water equation to obtain the tidal field; b) solves the gentle slope equation to calculate the wave field; c) realizes bidirectional wave-current coupling through the radiation stress gradient to generate a joint dataset containing the flow velocity vector, wave height field, and radiation stress field.

[0077] S7: Configure wave incident boundary conditions based on wave directional spectrum data, synthesize tidal boundary conditions based on tidal harmonic constants, and combine tidal levels with design wave heights to generate tidal wave combination conditions. Specifically, boundary conditions are constructed by: a) configuring incident boundary spectrum energy distribution based on wave directional spectrum data in S1; b) synthesizing tidal boundary time history curves based on astronomical tidal harmonic constants; c) combining extreme tidal levels with design wave heights to generate multiple conditions.

[0078] S8. Calculate the pile foundation stress parameters based on the wave-current combined field dataset generated in S6 and the tidal wave combined working condition in S7. Specifically, the combined dataset in S6 is input into the improved Morison equation: a) the resistance coefficient is dynamically adjusted according to the flow direction-pile angle (linearly increasing from 30° to 90°); b) the inertia force coefficient is corrected based on the turbulent kinetic energy intensity, and the shear force, bending moment extreme values, and fatigue spectrum of each cross section of the pile foundation are output;

[0079] S9, based on the wave-current combined field dataset generated by S6 and the tidal wave combined conditions of S7, uses the resonance characteristic equation to calculate the frequency matching and output the resonance amplification factor. Specifically, the harbor resonance assessment: a) Establish the harbor three-dimensional characteristic equation to solve the natural frequency; b) Analyze the matching degree between the main frequency of the wave spectrum in S6 and the natural frequency; c) Calculate the resonance amplification factor (range 0.1-3.0) based on the frequency overlap;

[0080] S10. Generate a structural safety assessment report based on pile foundation stress parameters and resonance amplification factors. Specifically, the pile foundation mechanical parameters from S8 and the resonance amplification factors from S9 are combined to generate a structural safety report using the Analytic Hierarchy Process (AHP) that includes the following elements: the ratio of maximum pile foundation stress to material yield strength; fatigue damage accumulation index; harbor resonance risk level; and a prioritized list of recommended reinforcement measures.

[0081] The above-mentioned technical solution comprehensively acquires a variety of key data, including seafloor topography, shoreline point clouds, wind speed and direction, current profiles, and wave directional spectrum data, providing a rich and comprehensive information foundation for subsequent analysis. Spatial registration is used to generate a seafloor digital elevation model, as well as time-series data processing and anomaly correction. This improves data quality and usability, ensuring greater accuracy in subsequent models and calculations based on this data. The construction of an unstructured grid and nested spatial discretization model, along with the rational setting of wave directional spectrum discretization parameters, enables more accurate simulation of complex hydrodynamic environments. The use of a variational assimilation algorithm optimizes the combination of hydrodynamic parameters, improving the fit to actual hydrodynamic conditions. Solving the wave-current coupling equation and calculating the radiation stress gradient correction term more realistically reflects the combined effects of waves and currents. The construction of combined tidal and wave conditions provides diverse and realistic scenarios for analyzing the stresses and other factors of port projects under different conditions. Calculating pile foundation stress parameters and resonance amplification factors, and generating a structural safety assessment report, provides a comprehensive and accurate assessment of port project safety, assisting in project decision-making and optimization, reducing project risks, and improving the quality and reliability of project construction.

[0082] Specifically, multi-source data collection and fusion are available. The shoreline point cloud data collection density can be set to 5 points / ㎡, 10 points / ㎡, or 15 points / ㎡. Vertical stratification of multi-layer wind speed data can be set to 10m, 20m, and 50m height layers. Wave directional spectrum sampling frequencies can be configured to 1Hz, 2Hz, or 5Hz. Seabed topography data can be collected using a multibeam echosounder (such as the Kongsberg EM series) or a side-scan sonar (such as the EdgeTech 4125). Shoreline point cloud collection can be performed using a terrestrial 3D laser scanner (such as the FARO Focus) or a drone-mounted LiDAR (such as the DJI L1). Current profile measurements can be performed using an ADCP (such as the Teledyne RDI Workhorse) or an electromagnetic current meter (such as the Nortek Aquadopp). Data storage media can be industrial-grade solid-state drives (such as the Samsung 870QVO). Sensor housings can be made of 316L stainless steel or IP68-rated waterproof materials. The multibeam echosounder was installed on the keel of the survey vessel, the shoreline LiDAR device was mounted beneath the drone's gimbal, and the ADCP was fixed to the top of the seafloor observation platform. Spatial registration used the ICP algorithm to align the point cloud data with the terrain data coordinate system. Time alignment was achieved through GPS clock synchronization, with an error limit of ±10ms. Anomaly correction used a Savitzky-Golay filter with a window length of 15 sampling points. Data fusion accuracy was improved to 0.1m, reducing temporal alignment error by 85%.

[0083] Nested dynamic mesh modeling and parameter optimization are implemented. The terrain gradient thresholds for the core infill area can be set to 0.15, 0.2, or 0.25, and the virtual cell interpolation weight coefficients can be set to 0.3, 0.5, or 0.7. The dynamic infill trigger threshold is set to a velocity gradient greater than 0.1s-1 or a sudden change in wave height greater than 0.5m. ANSYS Meshing or COMSOL Multiphysics can be used for mesh generation, and an NVIDIA A100 GPU cluster can be used for parameter optimization. Terrain gradient analysis is performed using the ArcGIS Pro spatial analysis module. The mesh model is stored in HDF5 or NetCDF format, and infill area marker data is encoded using binary bitmaps. Global mesh computing nodes are deployed on cloud computing servers, local infill modules are embedded in FPGA accelerator cards, and terrain gradient analysis components are integrated into GIS workstations. Terrain gradient thresholds are determined through historical data regression analysis, virtual cell interpolation weights are derived using the Kriging spatial interpolation method, and dynamic infill trigger conditions are set based on the Kolomogorov scaling theory. This improves grid computing efficiency by 40% and accelerates parameter optimization convergence by three times.

[0084] For wave-current coupling and engineering safety assessment, the Morison equation time step can be set to 0.05s, 0.1s, and 0.2s, the resonance frequency matching threshold can be set to 0.02Hz, 0.05Hz, and 0.1Hz, and the safety level classification threshold can be set to 0.6, 0.8, and 1.0. The wave-current coupling solver uses MIKE 21SW or TELEMAC, structural mechanics analysis uses ANSYS Mechanical or ABAQUS, and resonance characteristics are solved using the COMSOL Acoustics Module. A parameter library for pile foundation materials, such as Q345 steel, C50 concrete, and GFRP composites, is available. Viscous dampers (Taylor Devices Inc.) or tuned mass dampers can be used for damper selection. A loose coupling strategy is used for wave-current coupling iteration, with synchronization every five wave cycles. Resonant mode solutions use the Lanczos algorithm with a cutoff frequency of 1.0Hz. Safety assessment report generation is automatically populated based on an SQL database template. The calculation error of pile foundation force is controlled within ±5%, and the accuracy of resonance risk warning is increased to 92%.

[0085] In another technical solution, step S4 specifically includes:

[0086] S4.1. Extract the terrain gradient field based on the seafloor digital elevation model and divide the grid density control area according to the gradient threshold. Define the core density area (gradient ≥ 0.15), transition area (0.05 ≤ gradient < 0.15), and background sparse area (gradient < 0.05).

[0087] S4.2. Use the Delaunay triangulation algorithm to generate the initial unstructured grid, and superimpose the frontier advancing method in the core encryption area to achieve local secondary encryption, so that the grid side length meets the following requirements: core area ≤ 2m, transition area ≤ 5m, background area ≤ 10m;

[0088] S4.3. Construct a three-level nested hierarchy: a first-level global grid covering the entire computational domain with a resolution of 10m; a second-level local grid covering the port waters with a resolution of 5m; and a third-level refined grid focusing on the area surrounding the pile foundation cluster with a resolution of 2m.

[0089] S4.4. Couple nested boundaries via bidirectional flux interpolation and use virtual cells to transfer gradient fluxes of hydrodynamic parameters, including flow velocity and wave height;

[0090] S4.5. Based on the hydrodynamic parameter combination optimized in S5, dynamic adaptive encryption is triggered in the velocity shear layer (shear rate > 0.1s-1) and the wave breaking zone (wave steepness > 0.03), and the grid topology is adjusted in real time.

[0091] In this technical solution, by dividing grid density control zones based on terrain gradients, different grid densities can be set for varying degrees of terrain complexity. Core density is focused on areas of dramatic terrain change, while computational resources are rationally allocated to transition zones and sparse background areas, improving overall computational efficiency and accuracy. Delaunay triangulation combined with the frontier advancing method generates high-quality unstructured meshes that meet the mesh edge length requirements of different regions, ensuring the model's adaptability to complex terrain. A three-level nested hierarchy, with different levels tailored to different regional resolution requirements, ensures comprehensive coverage of the computational domain while achieving high-precision simulation in critical areas, effectively balancing computational effort and accuracy. A bidirectional flux interpolation method, coupled with nested boundaries and virtual cells to transfer hydrodynamic parameter gradient fluxes, ensures accurate transfer and sharing of hydrodynamic parameters between grid levels, improving overall model consistency and accuracy. Dynamic adaptive density in the velocity shear layer and wave breaking zone allows real-time adjustment of the grid topology based on actual hydrodynamic conditions, more accurately capturing complex hydrodynamic phenomena and significantly improving the accuracy of simulations of combined tidal and wave interactions, providing more reliable data support for port engineering construction.

[0092] Specifically, when implementing step S4, first, based on the generated seabed digital elevation model, a geographic information system (GIS) or professional terrain analysis software is used to extract the terrain gradient field. These software usually have powerful spatial analysis functions and can calculate the terrain gradient of each location based on the digital elevation model. For example, by calculating the elevation difference between adjacent grid points and combining the distance information, the terrain gradient value can be obtained. According to the set gradient threshold, the entire area is divided into a core density area (gradient ≥ 0.15), a transition area (0.05 ≤ gradient < 0.15) and a background sparse area (gradient < 0.05). In the actual port project, the reef area close to the coast is divided into a core density area due to the large terrain undulations. The terrain gradient is calculated to be greater than 0.15, and the area far away from the coast and with relatively flat terrain has a gradient less than 0.05, becoming a background sparse area.

[0093] Then, mesh generation is performed. The Delaunay triangulation algorithm is used to generate the initial unstructured mesh. This algorithm has good geometric properties, and the generated triangular mesh can well adapt to complex terrain boundaries. After the initial mesh is generated, the frontier advancing method is used for local secondary refinement in the core refinement area. The frontier advancing method starts from the existing mesh boundary and gradually advances inward to generate new mesh units. During the advancement process, the size and shape of the newly generated mesh are continuously adjusted according to the requirement that the mesh side length in the core refinement area must be ≤2m. For the transition zone, the mesh side length is controlled to ≤5m, and the mesh side length in the background zone is controlled to ≤10m. In actual operation, the mesh generation process is controlled through programming, and the corresponding parameters are set to ensure that the mesh side length meets the requirements.

[0094] When constructing the three-level nested hierarchy, the first level is the global grid, which covers the entire computational domain and has a resolution of 10m. This level of grid provides a macroscopic framework for the entire area and serves as a foundation for subsequent local refinement. The second-level local grid covers the port waters and has a resolution of 5m. Port waters are typically areas with complex tidal and wave effects. Increasing the grid resolution in this area helps to more accurately simulate hydrodynamic phenomena. For example, in the port's inlet and outlet channels, ships frequently sail, which has a significant impact on water flow and waves. The higher resolution of the second-level grid can better capture these changes. The third-level refined grid focuses on the area surrounding the pile group, with a resolution of 2m. The water flow and waves around the pile group are blocked and disturbed by the pile foundations, resulting in a very complex situation that requires extremely high grid resolution for accurate simulation. By rationally setting the resolution of the grids at different levels, we can ensure coverage of the entire computational domain while achieving high-precision simulation in key areas.

[0095] In order to realize nested boundary coupling, a bidirectional flux interpolation method is adopted. This method sets virtual cells at the boundaries of different levels of grids and uses these virtual cells to transfer the gradient flux of hydrodynamic parameters such as flow velocity and wave height. In the specific implementation process, according to the node positions and parameter values of the adjacent level grids, the parameter values at the virtual cells are obtained by interpolation calculation, thereby realizing the transfer of parameters between different levels of grids. For example, at the boundary between the first-level global grid and the second-level local grid, the flow velocity and wave height values near the boundary nodes in the first-level grid and the second-level grid are linearly interpolated to obtain the flow velocity and wave height at the virtual cells, and then these values are transferred to the adjacent grid levels, ensuring the continuity and accuracy of the hydrodynamic parameters between different levels of grids.

[0096] Finally, according to the hydrodynamic parameter combination optimized by S5, dynamic adaptive encryption is triggered in the velocity shear layer (shear rate>0.1s-1) and the wave breaking zone (wave steepness>0.03). By real-time monitoring of the hydrodynamic parameters, when it is found that the shear rate in a certain area is greater than 0.1s-1 or the wave steepness is greater than 0.03, the dynamic adaptive encryption program is started. The program increases the number of grids in the corresponding area and adjusts the grid topology according to pre-set rules. For example, in the velocity shear layer area, new grid lines may be inserted in the shear layer direction on the basis of the original grid to make the grid denser, thereby more accurately simulating the changes in flow velocity. Through this dynamic adaptive encryption method, the grid can be adjusted in real time according to the actual hydrodynamic conditions to improve the accuracy of the simulation.

[0097] In another technical solution, the step S4 specifically further includes:

[0098] S4.6. Extract the shoreline main axis direction θ based on the seabed digital elevation model in S2, use the range of θ ± 30° as the wave incident main direction interval, discretize the direction spectrum at 5° intervals within the wave incident main direction interval, and discretize at 15° intervals in the rest of the area;

[0099] S4.7. Based on the wave directional spectrum energy distribution collected in S1, three frequency bands are divided:

[0100] Core frequency band (energy share ≥ 70%): frequency interval Δf1 = 0.02 Hz;

[0101] Transition band (20% ≤ energy share < 70%): Δf2 = 0.05 Hz;

[0102] Edge frequency band (energy share <20%): Δf3 = 0.1 Hz;

[0103] S4.8. Hierarchical configuration of discrete parameters in three-level nested grids:

[0104] Level 1 global grid: direction interval 15°, frequency interval Δf3;

[0105] Secondary local grid: direction interval 10°, frequency interval Δf2;

[0106] Level 3 encrypted grid: direction interval 5°, frequency interval Δf1;

[0107] S4.9. Introduce a wavenumber-water depth coupling constraint. When the local water depth h and wavelength L of the DEM satisfy h / L ≤ 0.05, trigger the diffraction enhancement mode of the directional spectrum and increase the directional resolution of the area to 2°.

[0108] S4.10. Based on the velocity field assimilated in S5, activate the multi-directional wave interference compensation algorithm for areas where the velocity-wave angle is greater than 60°, and superimpose the reverse secondary wave component (with an intensity of 15%-30% of the forward wave) on the original directional spectrum.

[0109] In the above technical solution, by dividing the main wave incident direction interval according to the direction of the main axis of the shoreline and using discrete directional spectra with different intervals, we can focus on the main direction area of wave incidence, capture the wave direction information more carefully, and improve the simulation accuracy of the main wave propagation direction. For other areas, relatively coarse discrete intervals are used to reasonably balance the amount of calculation. According to the energy distribution of the wave directional spectrum, the three-level frequency band is divided and different frequency intervals are set. The core frequency band with a high energy proportion can be simulated in detail, and a wider frequency interval is used for the edge frequency band with lower energy, which effectively utilizes computing resources and improves simulation efficiency and accuracy. Discrete parameters are configured hierarchically in the three-level nested grid to adapt to the resolution and importance of grids at different levels, so that the model can simulate wave characteristics more accurately both globally and locally. The introduction of wave number-water depth coupling constraint to trigger the diffraction enhancement mode can improve the directional resolution under specific water depth-wavelength conditions, better simulate the wave diffraction phenomenon, and enhance the model's simulation ability of wave propagation in complex terrain. For areas where the velocity-wave angle is greater than 60°, the multi-directional wave interference compensation algorithm is activated and the reverse secondary wave component is superimposed. This takes into account the complex flow-wave interaction, can more realistically simulate the wave interference phenomenon, improve the accuracy of wave simulation in complex hydrodynamic environments, and provide more practical simulation results for port engineering construction.

[0110] Specifically, when implementing step S4, the operation of S4.6 is performed first. Based on the generated seabed digital elevation model (DEM), professional geospatial analysis software, such as ArcGIS, is used. The geometric shape of the coastline is analyzed through the vector analysis tool in the software to extract the direction θ of the main axis of the coastline. For example, multiple coordinate points of the coastline are fitted to obtain an axis representing the main trend of the coastline, thereby determining its direction. The range of θ±30° is used as the main direction interval of wave incidence. This is because in the actual marine environment, waves often have a major impact on port projects within this certain angle range related to the coastline. Within the main direction interval of wave incidence, a numerical calculation program is used to discretize the directional spectrum at intervals of 5°. In other areas outside the interval, it is discretized at intervals of 15°. This process is implemented through programming, and each discrete direction is numbered and recorded for subsequent use in the model.

[0111] When performing S4.7, the wave directional spectrum energy distribution data collected by S1 is analyzed. First, the proportion of wave energy in different frequency ranges is counted. For example, through spectrum analysis methods such as Fourier transform, the wave directional spectrum data is converted to the frequency domain, and then the sum of the energy in each frequency interval is calculated and compared with the total energy to obtain the energy proportion. For frequency ranges with an energy proportion ≥ 70%, they are defined as core frequency bands, and their frequency interval Δf1 is set to 0.02Hz. For energy proportions between 20% and 70%, they are defined as transition bands, and the frequency interval Δf2 is 0.05Hz. Areas with an energy proportion of less than 20% are edge bands, and the frequency interval Δf3 is set to 0.1Hz. In actual operation, data processing software, such as Python's data analysis library Pandas and scientific computing library NumPy, are used to write programs to implement the division of frequency bands and the setting of frequency intervals.

[0112] In S4.8, the discrete parameters are configured hierarchically for the three-level nested grids that have been constructed. For the first-level global grid, due to its wide coverage, the accuracy requirements are relatively low, and the direction interval is set to 15° and the frequency interval is Δf3. When programming, in the code module that defines the grid parameters, the direction and frequency related parameters of the first-level grid are set to the corresponding values. The second-level local grid covers the port waters and has higher accuracy requirements. The direction interval is set to 10° and the frequency interval is set to Δf2. The corresponding parameters are also set through the code. The third-level encrypted grid focuses on the area around the pile foundation group and has extremely high accuracy requirements. The direction interval is set to 5° and the frequency interval is set to Δfi. Ensure that the grids at each level are discretely set according to the set parameters to meet the requirements of wave simulation accuracy in different areas.

[0113] When implementing S4.9, the wave number-water depth coupling constraint is introduced. The local water depth h of each grid cell in the DEM is calculated, and the wavelength L of the corresponding position is calculated according to the wave theory formula. For example, for deep water waves, the wavelength L can be calculated by the formula L = gT2 / (2π), where g is the acceleration of gravity and T is the wave period. For shallow water waves, there is a corresponding formula for calculating the wavelength of shallow water waves. When h / L≤0.05, it is determined that the diffraction enhancement condition is met. Using the numerical calculation program, when an area that meets the conditions is detected, the directional resolution of the area is automatically increased to 2°. By modifying the relevant parameters of the directional discretization in the model, the directional resolution can be adjusted to better simulate the diffraction phenomenon of waves in shallow water areas or complex terrain.

[0114] In S4.10, the angle between the velocity and the waves is analyzed based on the velocity field data assimilated from S5. Computational fluid dynamics software is used to calculate the angle between the velocity vector and the wave propagation direction vector at each grid point through vector operations. When the angle is greater than 60°, the multi-directional wave interference compensation algorithm is activated. A reverse secondary wave component is superimposed on the original directional spectrum, with an intensity between 15% and 30% of the forward wave. The specific intensity is determined based on actual conditions by setting a random number generator in the program to randomly generate a value within the range of 15% to 30% as the intensity ratio of the reverse secondary wave component to the forward wave. In this way, complex flow-wave interactions are taken into account, enabling the model to more realistically simulate wave interference phenomena.

[0115] In another technical solution, step S5 specifically includes:

[0116] S5.1. Constructing a two-stage variational assimilation framework:

[0117] Tidal assimilation stage: Using the velocity profile data of S3 as a constraint, the bottom friction coefficient (range 0.001-0.003) and vertical eddy viscosity coefficient (range 0.1-1.2m 2 / s);

[0118] Wave assimilation stage: Using wave directional spectrum data as constraints, optimize the JONSWAP spectrum peak enhancement factor γ (range 1.5 to 3.3) and the breaking index κ (range 0.6 to 0.8);

[0119] S5.2. Design a terrain correction cost function, which includes the weighted sum of the following three terms:

[0120] The mean square error term between the observed flow rate and the flow rate calculated by the model;

[0121] The coupling constraint term between seafloor topography gradient and bottom shear stress;

[0122] Parameter sparsity regularization term;

[0123] The weight coefficient of the terrain gradient term is 0.3 to 0.6 times the weight of the observation term;

[0124] S5.3. Use the quasi-Newton-L-BFGS algorithm to solve the cost function. When the gradient descent rate is less than 5% for three consecutive iterations, the step size is automatically reduced to 1 / 4 of the original value.

[0125] S5.4. For areas where velocity and wave data conflict (relative error > 30%), dynamic weighting is applied according to the following rules:

[0126] Time dimension: centered on the moment of conflict initiation, with weight decaying exponentially with a half-life of 2 hours;

[0127] Spatial dimension: The weight of the area outside the 50-meter radius of the pile is reduced to 25% of the central area;

[0128] S5.5. Output the optimized parameter combination and confidence interval. When the coefficient of variation of tidal parameters is greater than 15% or the coefficient of variation of wave parameters is greater than 25%, manual review is triggered.

[0129] In this technical solution, a two-stage variational assimilation framework is constructed to optimize parameters for tides and waves separately. This fully utilizes different types of data, enabling the model to more accurately reflect the characteristics of tides and waves, as well as their interactions. A multi-factor terrain correction cost function is designed, with appropriate weighting coefficients. This comprehensively considers the accuracy of observational data, the impact of terrain on hydrodynamics, and the rationality of parameters, improving the reliability and adaptability of the model. A quasi-Newton-L-BFGS algorithm is used to solve the cost function, with the step size dynamically adjusted based on the gradient descent rate. This improves solution efficiency, accelerates convergence, and ensures optimized results within a reasonable timeframe. Dynamic weighting is applied to conflicting regions, with appropriate attenuation in both temporal and spatial dimensions, minimizing the negative impact of conflicting data on the assimilation results and ensuring they are more realistic. The optimized parameter combinations and confidence intervals are output, and a manual review mechanism is implemented to ensure parameter accuracy and reliability. When the coefficient of variation of a parameter exceeds a reasonable range, prompt manual intervention is implemented to prevent erroneous results from misleading subsequent engineering analysis, providing a more reliable basis for hydrodynamic parameters in port engineering construction.

[0130] Specifically, when implementing step S5, first perform the operation of S5.1. Construct a two-stage variational assimilation framework. In the tidal assimilation stage, use numerical simulation software such as COMSOL Multiphysics or FVCOM (finite volume coastal ocean model). Import the processed velocity profile data obtained in S3 into the software. These data contain velocity information at different depths, different positions and different times. Set the optimization target in the software, that is, by adjusting the bottom friction coefficient (the value range is between 0.001-0.003) and the vertical eddy viscosity coefficient (the value range is between 0.1-1.2m 2 / s), so that the velocity calculated by the model is as close as possible to the actual observed velocity profile data. For example, in a tidal simulation scenario, the initial bottom friction coefficient is set to 0.002 and the vertical eddy viscosity coefficient is set to 0.5m 2 / s, through continuous iterative calculation, these two parameters are adjusted according to the feedback of velocity profile data to reduce the error between the model calculated velocity and the observed velocity.

[0131] In the wave assimilation stage, the above-mentioned numerical simulation software is also used. The collected wave directional spectrum data is input into the software and used as a constraint to optimize the JONSWAP spectrum peak enhancement factor γ (with a value range of 1.5-3.3) and the breaking index κ (with a value range of 0.6-0.8). The wave directional spectrum data contains wave energy distribution information of different directions and frequencies. By adjusting the values of γ and κ, the wave characteristics simulated by the model, such as wave height and period, match the actual observed wave directional spectrum data. For example, when simulating the wave conditions in a sea area, the initial setting of γ is 2.0 and κ is 0.7. Based on the comparative analysis of the wave directional spectrum data, these two parameters are continuously adjusted to improve the model's simulation accuracy of the waves.

[0132] When performing S5.2, design a terrain correction cost function. Within the programming environment of the numerical simulation software, write code to construct the cost function. The cost function consists of three weighted sums. The first term is the mean squared error (MSE) between the observed and model-calculated velocities. This term is calculated by summing the squares of the differences between the observed and model-calculated velocities at each observation point and averaging them. The second term is a coupled constraint between the seafloor topography gradient and the bottom shear stress. First, terrain gradient information is extracted from the seafloor digital elevation model using Geographic Information System (GIS) software. Then, based on the principles of fluid mechanics, the bottom shear stress is calculated. This coupled constraint is then used to calculate the value of this term. The third term is a parameter sparsity regularization term, which ensures that the optimized parameters are not overly complex and have a certain degree of sparsity. When setting the weight coefficients, set the terrain gradient term to 0.3-0.6 times the observation term weight. For example, if the observation term weight is set to 1, the terrain gradient term weight can be set to 0.4. By adjusting the weight coefficients, the contribution of each term in the cost function is balanced.

[0133] In S5.3, the quasi-Newton-L-BFGS algorithm is used to solve the cost function. The numerical simulation software uses a relevant optimization algorithm library, such as the optimization function in the SciPy library in Python. The constructed terrain correction cost function is used as input, and the algorithm's initial parameters, such as the initial step size, are set. During the iterative calculation process, the gradient descent rate is monitored in real time. When the gradient descent rate for three consecutive iterations is less than 5%, the step size is automatically reduced to 1 / 4 of the original value through programming. For example, the initial step size is set to 0.1. When the gradient descent rate condition is met, the step size is adjusted to 0.025. The iterative calculation continues until the convergence condition is met, resulting in the optimized parameter value.

[0134] When implementing S5.4, the conflicting areas of velocity and wave data are processed. In the numerical simulation software, the velocity and wave related data calculated by the model are compared by programming. When the relative error is greater than 30%, it is determined to be a data conflict area. In the time dimension, the weight is set using an exponential decay function with the conflict starting time as the center. For example, assuming the conflict starting time is t0 and the current time is t, the weight The half-life is 2 hours. In the spatial dimension, the distance from each point to the pile is calculated, with the pile as the center. When the distance is greater than 50 meters, the weight is reduced to 25% of the central area. This reduces the impact of conflicting data during the data assimilation process, making the assimilation results more reliable.

[0135] In S5.5, after parameter optimization is completed, an output function is set in the numerical simulation software to output the optimized parameter combination and confidence interval. For example, the output parameters such as the bottom friction coefficient of 0.0025±0.0003 and the vertical eddy viscosity coefficient of 0.6±0.05 and their confidence intervals are output. At the same time, a judgment condition is set. When the coefficient of variation of the tidal parameters is greater than 15% or the coefficient of variation of the wave parameters is greater than 25%, the manual review mechanism is triggered. The corresponding judgment code is written in the software. When it is detected that the parameter variation coefficient is out of range, a prompt message is issued to remind researchers to conduct manual inspection and analysis to ensure the accuracy of the parameters.

[0136] In another technical solution, step S6 specifically includes:

[0137] S6.1. Solve the wave-current coupling equations in stages:

[0138] Wave-dominated stage: Based on the discrete directional spectrum parameters of S4, the gentle slope equation is used to calculate the wave field and output the wave height, wave direction and breaking energy distribution;

[0139] Tide-dominated stage: Based on the hydrodynamic parameters optimized by S5, the three-dimensional shallow water equations are solved to obtain the flow velocity field and water level field;

[0140] Joint correction stage: the wave radiation stress gradient is injected into the tidal current equation as a source term, and the velocity field is fed back into the wave refraction calculation;

[0141] S6.2. Calculate the radiation stress gradient correction term:

[0142] In the wave breaking zone (wave steepness > 0.03), a turbulent kinetic energy compensation model is used to convert 10%-15% of the breaking dissipated energy into radiation stress increments;

[0143] The step length is calculated by dynamically adjusting the correction term based on the terrain gradient of S2. The step length in steep slope areas (gradient ≥ 0.1) is ≤ 2m, and the step length in flat areas is ≤ 10m.

[0144] S6.3. Perform bidirectional coupling iteration:

[0145] After each wave-current coupling calculation is completed, the relative error of wave height (threshold ≤ 3%) and flow velocity deviation (≤ 0.05 m / s) between two adjacent calculations are compared;

[0146] The calculation is terminated when the threshold is not exceeded for three consecutive iterations, otherwise the radiation stress field is updated and iterated again;

[0147] S6.4. Generate a wave-current joint field dataset, including:

[0148] Space vector field: flow velocity, flow direction, wave height, wave direction;

[0149] Energy scalar field: radiation stress intensity, turbulent kinetic energy density;

[0150] Dynamic characteristic field: location of the broken zone and vortex core area.

[0151] In this technical solution, the wave-current coupling equations are solved in stages, making the calculation of wave and tidal fields more targeted. The appropriate equations and parameters for each are utilized to improve computational accuracy. The wave-dominated stage accurately simulates wave characteristics, while the tidal stage accurately obtains flow velocity and water level information. The joint correction stage simulates the interaction between the two, fully reflecting the wave-current coupling phenomenon. When calculating the radiation stress gradient correction term, a reasonable model is used to convert energy in the wave breaking zone, and the step size is dynamically adjusted based on the terrain gradient. This more accurately considers the influence of wave breaking and topographic factors on waves and currents, improving the realism of the simulation. Bidirectional coupling iterations are performed, strictly controlling the relative error in wave height and flow velocity deviation to ensure stable and accurate calculation results and avoid deviations due to error accumulation. This generates a combined wave-current field dataset containing comprehensive information, providing a rich, comprehensive, and accurate data foundation for subsequent structural stress analysis of port projects, assisting project design and construction and improving the port project's ability to cope with complex hydrodynamic environments.

[0152] Specifically, when implementing step S6, the work of S6.1 is carried out first. In the wave-dominated stage, numerical simulation software, such as the SWAN (Simulating WAves Nearshore) model, is used. The discrete directional spectrum parameters determined in S4 are input into the SWAN model. These parameters include information such as the wave energy distribution in different directions and frequencies. The wave field is calculated using the gentle slope equation. The gentle slope equation can take into account phenomena such as refraction and diffraction of waves during propagation, and more accurately simulate the propagation of waves on complex terrain. Through model calculation, wave height, wave direction and breaking energy distribution data are output. For example, in a sea area simulation with complex seabed topography, the wave height at different locations is calculated by the SWAN model. The wave height in some areas can reach 3 meters, and the wave direction shows obvious changes affected by the terrain. At the same time, the energy distribution in the wave breaking area is determined.

[0153] In the tidal-dominated stage, numerical models such as FVCOM (Finite-Volume Community Ocean Model) are used. The hydrodynamic parameters obtained by S5 optimization, such as the bottom friction coefficient and the vertical eddy viscosity coefficient, are input into the FVCOM model. The three-dimensional shallow water equation is used to solve the problem. The three-dimensional shallow water equation takes into account factors such as the movement of water in the horizontal and vertical directions and the change of water depth. Through model calculation, the flow velocity field and water level field information are obtained. In the tidal simulation of a certain port, the FVCOM model calculates the flow velocity at different positions and depths. The flow velocity in the center of the channel can reach 1.5m / s. At the same time, the water level changes in the entire port area are obtained.

[0154] During the joint correction phase, the wave radiation stress gradient is injected as a source term into the tidal current equation, implemented through programming within the numerical simulation software. Simultaneously, the velocity field is fed back into the wave refraction calculation to adjust parameters such as the wave propagation direction. For example, in a simulation of an estuary, the influence of the wave radiation stress gradient on the tidal current was taken into account, resulting in significant changes in tidal velocity and direction. Simultaneously, the changes in the velocity field also resulted in different wave refraction, resulting in a more realistic simulation of wave-current interaction.

[0155] When performing S6.2, a turbulent kinetic energy compensation model is used in wave breaking areas (wave steepness > 0.03). This model is implemented using computational fluid dynamics software, such as OpenFOAM. 10%-15% of the dissipated energy from breaking is converted into radiation stress increments, and an appropriate value, such as 12%, is selected within this range based on the actual situation. The correction term calculation step size is dynamically adjusted based on the terrain gradient of S2, and terrain gradient information is extracted from the seabed digital elevation model using geographic information system (GIS) software. In steep slope areas (gradient ≥ 0.1), the calculation step size is set to no more than 2 meters through programming, and to no more than 10 meters in gentle areas. For example, in an area with large seabed topography, the step size is set to 1.5 meters on steep slopes and 8 meters in gentle areas to ensure calculation accuracy.

[0156] In S6.3, a bidirectional coupling iteration is performed. An iterative calculation program is set up in the numerical simulation software. After each wave-current coupling calculation is completed, the relative error of wave height and the flow velocity deviation of two adjacent times are automatically compared. When the wave height relative error threshold is ≤3% and the flow velocity deviation is ≤0.05m / s, the next iteration is continued. When these thresholds are not exceeded for three consecutive iterations, the calculation results are judged to be converged and the calculation is terminated. If the conditions are not met, the radiation stress field is updated and the iterative calculation is performed again. For example, in a certain simulation, the relative errors of wave height in the first two iterations were 4% and 3.5% respectively, and the flow velocity deviations were 0.06m / s and 0.055m / s respectively. The conditions were not met. After updating the radiation stress field, the iteration continued. After multiple iterations, the termination conditions were finally met.

[0157] In S6.4, after completing the wave-current coupling calculation, the numerical simulation software sets the output function to generate a combined wave-current field dataset. The calculated data, such as flow velocity, flow direction, wave height, and wave orientation, are organized into a spatial vector field. Data such as radiation stress intensity and turbulent kinetic energy density are calculated to form an energy scalar field. Information such as the location of the fracture zone and the vortex core region is determined to construct a dynamic characteristic field. These different types of data are integrated to generate a complete combined wave-current field dataset, providing data support for subsequent analysis.

[0158] In another technical solution, step S7 specifically includes:

[0159] S7.1. Wave incident boundary configuration:

[0160] According to the measured wave direction spectrum of S1, the main direction θ and the spectrum peak period Tp are extracted, and the JONSWAP spectrum type incident wave is generated at the boundary of the calculation domain;

[0161] Based on the seabed topography gradient of S2, a wave refraction correction angle Δθ (range 1°~5°) is automatically added in nearshore shallow waters (water depth h<10m);

[0162] S7.2. Tidal boundary synthesis:

[0163] Eight main tides (M2, S2, K1, O1, P1, Q1, N2, and K2) are selected, and the tide level time series is generated based on the harmonic constants (amplitude H and lag angle g):

[0164] Full amplitude was retained for subtidal waves with >90% confidence;

[0165] The tidal amplitude with a confidence level of 60% to 90% is attenuated by 30% to 50%;

[0166] Tidal sub-segments with a confidence level <60% were eliminated;

[0167] A topographic attenuation factor was introduced to apply a 10% to 20% attenuation to the tidal amplitude in the closed area of the bay top (shoreline curvature > 0.1 / m);

[0168] S7.3. Constructing tide and wave combined working conditions:

[0169] Conventional operating conditions: 9 benchmark operating conditions are generated by orthogonal combination of tide type (spring tide, moderate tide, neap tide) and wave return period (1 year, 50 years);

[0170] Extreme conditions: Using the asynchronous extreme value method, we randomly select combinations of wave height (100-year return period), tide level (astronomical high tide + storm surge), and wind direction (opposite to the tidal current) to generate 50 random extreme scenarios.

[0171] S7.4, Boundary Phase Calibration:

[0172] When the time difference between the surge moment and the wave crest arrival time is greater than 15 minutes, a time offset is applied to the wave boundary to ensure the maximum wave-current coupling strength;

[0173] Within the breakwater shelter area (shielding angle ±30°), the wave incident direction is slightly adjusted by 2° to 3° to compensate for the diffraction effect.

[0174] In the above technical solution, the wave incident boundary configuration uses JONSWAP spectral incident waves generated based on the measured wave directional spectrum, which more realistically reflects actual wave characteristics. In nearshore shallow waters, a wave refraction correction angle is added based on the terrain gradient to accurately simulate the changes in wave propagation in complex terrain, improving wave simulation accuracy. When synthesizing tidal boundaries, the main tidal components are selected and the amplitudes are processed based on confidence and terrain characteristics to make the synthesized tidal time series more realistic. This takes into account the reliability of different tidal components and the impact of terrain on tidal levels. A wave and tide combination operating condition is constructed, with both conventional and extreme operating conditions, to comprehensively cover the different hydrodynamic conditions that ports may face, providing a foundation for multi-scenario analysis of port projects. Boundary phase calibration optimizes wave-current coupling through time offset and fine-tuning of the wave incident direction, ensuring that simulation results are closer to actual wave-current interaction. This provides more accurate hydrodynamic data for port project design and enhances the project's ability to cope with complex hydrodynamic environments.

[0175] Specifically, when implementing step S7, first perform the operation in step S7.1. Data analysis software, such as Python's pandas and numpy libraries, is used to process the wave directional spectrum data measured in step S1. Using a spectrum analysis algorithm, the principal direction θ and the spectral peak period Tp are extracted. For example, in a wave directional spectrum dataset containing multiple directions and frequencies, analysis determines that the principal direction is 30° northeast by east, and the spectral peak period is 8 seconds. In numerical simulation software, such as COMSOL Multiphysics or FVCOM (Finite Volume Coastal Ocean Model), the provided wave generation module is used to generate JONSWAP spectral incident waves at the computational domain boundary. This spectral type can effectively simulate the wave energy distribution in the actual ocean. Simultaneously, geographic information system (GIS) software is used to extract terrain gradient information from the seafloor digital elevation model in step S2. In nearshore shallow waters (water depth h < 10m), a wave refraction correction angle Δθ is automatically added through programming. An appropriate value is selected within the range of 1°-5° based on the actual terrain conditions. For example, in a shallow nearshore area with complex seabed topography, Δθ is set to 3° after analysis and calculation.

[0176] When operating S7.2, set up the tidal boundary synthesis module in the numerical simulation software. Eight major tidal components were selected: M2, S2, K1, O1, P1, Q1, N2, and K2. The harmonic constants (amplitude H and lag angle g) for these tidal components were obtained from previous ocean observation data. Using trigonometric relationships, tidal time series were generated based on the harmonic constants. For example, tidal level η can be calculated using the formula η = H·cos(ωt+g), where ω is the angular frequency of the tidal component, t is time, H represents the amplitude, and g represents the initial phase. For tidal components with a confidence level greater than 90%, the program is set to retain their full amplitude. For tidal components with a confidence level between 60% and 90%, code is written to attenuate their amplitude by 30% to 50%. Tidal components with a confidence level less than 60% are directly discarded during data processing. At the same time, GIS software was used to analyze shoreline curvature. For closed bay crest areas (shoreline curvature > 0.1 / m), a terrain attenuation factor was introduced into the tide calculation program, attenuating the tide amplitude by 10%-20%. For example, in a closed bay crest area, the tide amplitude was attenuated by 15%.

[0177] In S7.3, a tide and wave combination working condition is constructed. In the working condition setting module of the numerical simulation software, for conventional working conditions, orthogonal combinations are made according to the tide type (spring tide, moderate tide, neap tide) and the wave return period (1 year, 50 years). Through loop statements and conditional judgment statements, 9 benchmark working conditions are generated. For example, the combination of spring tide and 1-year wave return period, the combination of spring tide and 50-year wave return period, etc. For extreme working conditions, the asynchronous extreme value method is adopted. A random number generation function is written in the program to randomly extract combinations of wave height (once in a hundred years), tide level (astronomical spring tide + storm surge) and wind direction (opposite to the tide) to generate 50 random extreme scenarios. For example, a combination of a wave height of 10 meters (once in a hundred years), a tide level of astronomical spring tide plus 3 meters of storm surge, and a wind direction opposite to the tide is randomly extracted, and a sufficient number of extreme scenarios are generated through multiple random extractions.

[0178] When implementing S7.4, set up a boundary phase calibration module in the numerical simulation software. Monitor the surge moment and the wave crest arrival time in real time. When the deviation is greater than 15 minutes, apply a time offset to the wave boundary through programming. For example, if the surge moment is 20 minutes earlier than the wave crest arrival time, delay the wave boundary time by 20 minutes to ensure the maximum wave-current coupling intensity. Within the breakwater shelter area (shelter angle ±30°), use the direction adjustment function in the model to fine-tune the wave incident direction by 2°-3°. For example, within a breakwater shelter area, fine-tune the wave incident direction by 2.5° to compensate for the diffraction effect and make the simulation results more consistent with the actual situation.

[0179] In another technical solution, the pile foundation force calculation step S8 specifically includes:

[0180] S8.1. Input the wave-current combined field data set of S6 into the modified Morison equation at a time step of 0.1 seconds, and dynamically adjust the drag coefficient C. d (range 1.8~2.5) and inertia coefficient C m (Range 1.6~2.0):

[0181] When the angle between the flow direction and the pile axis is greater than 45°, C d Increase by 10% to 15%;

[0182] In the vortex shedding region (turbulent kinetic energy>0.3J / m 3 ), C m Reduce by 5% to 8%;

[0183] S8.2. Calculate the time history curves of shear force and bending moment for each section of the pile foundation and extract the maximum value, standard deviation, and fatigue damage index.

[0184] In the above technical solution, the wave-current combined field data set is input at a time step of 0.1 seconds, which can accurately capture the influence of the change of wave and current over time on the pile foundation stress and improve the calculation accuracy. Dynamic adjustment of the resistance coefficient C d and the inertia coefficient C m , taking into account the angle between the flow direction and the axis of the pile foundation and the influence of the vortex shedding zone, the calculation parameters are more in line with the actual interaction between waves and currents and the pile foundation, greatly improving the accuracy of the calculation results of the improved Morison equation. By calculating the shear force and bending moment time history curves of each section of the pile foundation and extracting the maximum value, standard deviation and fatigue damage index, the stress state of the pile foundation under the long-term action of waves and currents can be comprehensively evaluated. The maximum value can reflect the extreme load borne by the pile foundation, the standard deviation reflects the stress fluctuation, and the fatigue damage index is used to evaluate the durability of the pile foundation, providing comprehensive and critical data support for the design optimization, safety monitoring and maintenance of port engineering pile foundations, ensuring the long-term stability and safety of port engineering structures.

[0185] Specifically, when implementing step S8, first perform the operation in step S8.1. Import the wave-current combined field dataset generated in step S6 into professional structural mechanics analysis software, such as ANSYS or ABAQUS. These software programs possess powerful numerical computing capabilities and a rich library of physical models, enabling them to handle complex mechanical calculations. Set the calculation parameters in the software, including a time step of 0.1 seconds, to ensure that the dynamic changes of the wave-current combined field over time are captured in detail.

[0186] For the improved Morison equation, write the corresponding calculation module in the software or call the built-in related calculation function. This equation takes into account the drag force and inertia force of the water flow on the pile foundation, and its expression is generally Where F is the force per unit length of the pile, ρ is the fluid density, and C dis the resistance coefficient, A is the projected area of the pile foundation perpendicular to the flow direction, U is the flow velocity, C m is the coefficient of inertia, V is the volume of fluid displaced by the pile foundation, and dt / dU is the acceleration.

[0187] During the calculation process, the resistance coefficient C is dynamically adjusted according to the actual situation. d and the inertia coefficient C m . Using the conditional judgment statement of the software, when the angle between the flow direction and the pile axis is greater than 45°, the original C d The initial C d When the value is 2.0, when the angle between the flow direction and the pile axis meets the conditions, C d Adjusted to 2.0×(1+12%)=2.24. In the vortex shedding area, through the area recognition function of the software, when the turbulent kinetic energy is detected to be greater than 0.3J / m 3 When C m The value (range 1.6-2.0) decreased by 5%-8%. Assuming the initial C m is 1.8, and is adjusted to 1.8×(1-6%)=1.692 in the vortex shedding region.

[0188] When operating S8.2, the structural mechanics analysis software utilizes its structural analysis module to calculate the shear force and bending moment at each section of the pile foundation based on the forces calculated using the modified Morison equation. Numerical methods such as time integration generate time history curves of the shear force and bending moment at each section. The software automatically records the shear force and bending moment values for each section at each time step, generating a complete time history data set.

[0189] After obtaining the time history curve, the key indicators are extracted using the data processing function of the software. By searching for the maximum value in the time history data, the maximum value of the shear force and bending moment of each section of the pile foundation is obtained, which can reflect the maximum load borne by the pile foundation during the entire calculation process. The standard deviation of the time history data is calculated to evaluate the degree of fluctuation in the force on the pile foundation. The larger the standard deviation, the more severe the force fluctuation, and the greater the fatigue damage impact on the pile foundation may be. For the calculation of the fatigue damage index, Miner's linear cumulative damage theory is used. By statistically analyzing the stress cycles in the time history curve and combining the fatigue performance parameters of the pile foundation material, the fatigue damage index is calculated. For example, assuming that the SN curve of the pile foundation material is known, the fatigue damage index is calculated according to the Miner theory based on the stress amplitude and number of cycles in the time history curve. These key indicators are sorted out and output to provide a data basis for subsequent pile foundation design and safety assessment.

[0190] In another technical solution, step S9 specifically includes:

[0191] S9.1. Establish a three-dimensional characteristic equation based on the seabed digital elevation model (DEM) and the harbor contour line in S2, and solve for the first six resonant modes (frequency range 0.01-0.5 Hz);

[0192] S9.2. Compare the main frequency of the wave spectrum of S6 with the resonant mode frequency and calculate the matching degree:

[0193] When the main frequency difference is ≤0.02Hz, it is determined to be a strong resonance risk, and the amplification factor is 1.5 to 3.0;

[0194] When the main frequency difference is >0.1Hz, it is determined to be non-resonant and the amplification factor is set to 1.0;

[0195] S9.3. At the breakwater opening (width > 50 m), the effect of secondary reflected waves is superimposed, and an increment of 10% to 20% is added to the amplification factor.

[0196] In the above technical solution, a three-dimensional characteristic equation is established based on the DEM and the harbor contour line, and the resonance mode is solved. This can accurately analyze the frequency range of the harbor that may resonate under the action of waves from the perspective of terrain and harbor structure, providing basic data for subsequent risk assessment. By comparing the matching degree between the main frequency of the wave spectrum and the resonance mode frequency calculation, clarifying the judgment criteria for strong resonance risk and no resonance, and setting the corresponding amplification factor, the resonance risk level of the harbor under different wave conditions can be intuitively quantified, helping engineers to quickly judge the potential risks faced by the engineering structure. The influence of secondary reflected waves is superimposed at the breakwater opening and the amplification factor is incrementally corrected, fully considering the complex hydrodynamic phenomena in special areas of the port, making the resonance risk assessment results closer to the actual situation, and providing a more reliable basis for the design optimization, safety monitoring and maintenance decision-making of port projects, effectively improving the ability of port facilities to cope with complex marine environments, and ensuring the long-term stable operation of the port.

[0197] Specifically, when implementing step S9, first carry out the work of S9.1. Use professional engineering analysis software, such as COMSOL Multiphysics or ANSYS. Import the seabed digital elevation model (DEM) data generated by S2 into the software. These data record the elevation information of the seabed terrain in detail. At the same time, import the harbor contour line data. The harbor contour line clarifies the boundary shape and range of the harbor. In the software, based on the principles of fluid mechanics and structural dynamics, these data are used to establish a three-dimensional characteristic equation. Taking COMSOL Multiphysics as an example, in the model building module, select a suitable physical field interface, such as the acoustic-structural interaction interface, and set relevant parameters such as fluid density, elastic modulus, etc. according to the actual physical characteristics of the harbor, so as to establish a three-dimensional characteristic equation that accurately reflects the mechanical characteristics of the harbor under the action of waves.

[0198] After establishing the equations, the software's solver is used to solve for the first six resonant modes. In the solver settings, the frequency range is limited to 0.01-0.5 Hz, a common resonant frequency range determined based on practical engineering experience and research on harbor resonance phenomena. The solver iterates the three-dimensional characteristic equations using numerical calculation methods, such as the finite element method, to determine the frequencies of the first six resonant modes. For example, the first resonant mode is calculated to have a frequency of 0.05 Hz, the second at 0.12 Hz, and so on. These frequency values are then recorded for subsequent analysis.

[0199] When performing S9.2, extract the dominant frequency of the wave spectrum from the wave-current combined field data set of S6. In data analysis software, such as Python's numpy and scipy libraries, write code to perform spectral analysis on the wave spectrum data. Using algorithms such as Fourier transform, determine the frequency at which the wave spectrum has the most concentrated energy, i.e., the dominant frequency. Compare the extracted dominant frequency of the wave spectrum with the previously calculated resonant modal frequency. Write a comparison program to calculate the dominant frequency difference. When the dominant frequency difference is ≤0.02Hz, set the program to determine a strong resonance risk and select an amplification factor within the range of 1.5-3.0 based on the actual situation, such as 2.0. When the dominant frequency difference is >0.1Hz, determine that there is no resonance and set the amplification factor to 1.0. In this way, based on the relationship between the dominant frequency of the wave spectrum and the resonant modal frequency, the resonance risk status of the harbor can be quickly and accurately assessed.

[0200] In S9.3, Geographic Information System (GIS) software is used to analyze the location and opening width of the breakwater. When the width of the breakwater opening is detected to be greater than 50m, the impact of secondary reflection waves is considered in the resonance risk assessment procedure. Secondary reflection waves are formed by waves reflecting again after being reflected at the breakwater opening, which will enhance the effect of waves on the harbor. A 10%-20% increment is added to the amplification factor. For example, if the previously calculated amplification factor is 1.8, when the conditions are met at the breakwater opening, the amplification factor is adjusted to 1.8×(1+15%)=2.07. In this way, the amplification factor is corrected to make the resonance risk assessment results more consistent with the actual situation.

[0201] In another technical solution, the generation of the S10 security assessment report specifically includes:

[0202] S10.1. Construct a dual-index evaluation matrix:

[0203] Transverse dimension: ratio of the maximum bending moment of the pile foundation to the design bearing capacity (threshold 1.0);

[0204] Longitudinal dimension: product of resonance amplification factor and allowable vibration acceleration (threshold 0.8);

[0205] S10.2, divided into four levels of safety:

[0206] Level I (safe): Both indicators are <0.6, marked green;

[0207] Level II (warning): Any indicator is 0.6-0.8, marked yellow;

[0208] Level III (risk): Any indicator is 0.8 to 1.0, marked orange;

[0209] Level IV (dangerous): any indicator > 1.0, marked red;

[0210] S10.3. Automatically generate optimization suggestion library:

[0211] The orange level triggers an increase in pile diameter (recommended increase of 10% to 15%) or a damper layout plan;

[0212] The red level triggers measures such as modification of the harbor layout or resonant frequency shift.

[0213] In the above technical solution, a dual-index evaluation matrix is constructed, which comprehensively considers the relationship between the maximum bending moment of the pile foundation and the design bearing capacity, as well as the product of the resonance amplification factor and the allowable vibration acceleration. The safety of the engineering structure is comprehensively evaluated from two key aspects: the pile foundation force and the resonance effect, making the evaluation results more scientific and comprehensive. The four-level safety level is divided and equipped with intuitive color markings, which can clearly and intuitively show engineers and decision makers the safety status of the engineering structure, facilitating quick judgment and taking corresponding measures. An optimization suggestion library is automatically generated to provide targeted optimization suggestions according to different safety levels, such as increasing the pile diameter or laying out dampers at the orange level, and modifying the port basin plane layout or resonant frequency offset measures at the red level. This provides a clear direction for the optimization and improvement of the project, helps to improve the safety of the project, reduce potential risks, and ensure the stable operation and long-term development of the port project.

[0214] Specifically, when implementing step S10, first perform the operation of S10.1. Use data analysis software, such as Python's pandas library and numpy library, to construct a dual-index evaluation matrix. Obtain the maximum bending moment data and design bearing capacity data of the pile foundation from the previous calculation results, and calculate the ratio of the two as the lateral dimension data. For example, the maximum bending moment of a pile foundation calculated by structural mechanics analysis software is 800kN·m, and its design bearing capacity is 1000kN·m. The lateral dimension value of the pile foundation is 800÷1000=0.8. For the longitudinal dimension, obtain the resonance amplification factor and the allowable vibration acceleration data, and calculate the product of the two. Assume that the resonance amplification factor of a certain area is 1.5, and the allowable vibration acceleration is 0.5m / s 2 , the vertical dimension is 1.5 × 0.5 = 0.75. These calculated data are organized into a matrix for subsequent analysis.

[0215] When performing S10.2, a program is written in the data analysis software to classify the safety level into four levels based on the data in the dual-index evaluation matrix. When both indicators are less than 0.6, the program sets the marker to green, indicating Level I (safe). For example, if the lateral dimension value of a pile foundation is 0.5 and the longitudinal dimension value is 0.4, satisfying the dual-index value of <0.6, it is marked green. If either indicator is between 0.6 and 0.8, it is marked yellow, indicating Level II (warning). If the lateral dimension value is 0.7 and the longitudinal dimension value is 0.5, and the lateral dimension value is between 0.6 and 0.8, it is marked yellow. If either indicator is between 0.8 and 1.0, it is marked orange, indicating Level III (risk). If the longitudinal dimension value is 0.9 and the lateral dimension value is 0.7, and the longitudinal dimension value is between 0.8 and 1.0, it is marked orange. If either indicator is greater than 1.0, it is marked red, indicating Level IV (danger). In this way, different safety levels are clearly defined.

[0216] In S10.3, a program that automatically generates an optimization suggestion library is established using programming software, such as Python. When the safety level is detected to be orange, the program automatically generates suggestions for increasing the pile diameter (recommended increase of 10%-15%) or damper layout plans. For example, for a pile foundation at the orange level, the program generates suggestions: "The pile foundation is in level III (risk) status. It is recommended to increase the pile diameter with a recommended increase of 12%" or "It is recommended to lay dampers near the pile foundation to reduce the impact of vibration." When the safety level is detected to be red, the program generates suggestions for modifying the port basin plan layout or resonant frequency shift measures. For example, "This area is in level IV (dangerous) status. It is recommended to modify the port basin plan layout and adjust the port basin shape or position" or "It is recommended to take resonant frequency shift measures, such as changing the stiffness of the structure in the port basin to avoid the resonant frequency." Through these program settings, optimization suggestions can be automatically generated according to the safety level.

[0217] A tidal wave combined testing system for executing the method is provided, comprising:

[0218] A data acquisition module is used to obtain seabed topography data, shoreline point cloud data, multi-layer wind speed and direction data, flow velocity profile data, and wave direction spectrum data;

[0219] The elevation model generation module is used to spatially register terrain data with shoreline point cloud data to generate a seabed digital elevation model;

[0220] A data processing module is used to perform spatiotemporal alignment processing on the collected multi-layer wind speed and direction data, flow velocity profile data, and wave direction spectrum data to generate a quality verification data set with temporal and spatial consistency;

[0221] A model construction module is used to construct an unstructured grid based on the seabed digital elevation model, generate a nested spatial discrete model, and set wave directional spectrum discrete parameters in the nested spatial discrete model;

[0222] Parameter optimization module, which is used to optimize the S3 quality verification data set using the variational assimilation algorithm and output the optimized hydrodynamic parameter combination;

[0223] The wave-current combined field generation module is used to solve the wave-current coupling equation and calculate the radiation stress gradient correction term based on the combination of wave directional spectrum discrete parameters and hydrodynamic parameters to generate a wave-current combined field data set;

[0224] The working condition construction module is used to configure the wave incident boundary conditions according to the wave direction spectrum data, synthesize the tide level boundary conditions according to the tidal harmonic constant, and combine the tide level with the design wave height to generate the tidal wave combination working condition;

[0225] The pile foundation force calculation module is used to calculate the pile foundation force parameters based on the generated wave-current joint field data set and the tidal wave combined working condition;

[0226] The resonance factor calculation module is used to calculate the frequency matching degree using the resonance characteristic equation based on the generated wave-current joint field data set and the tide-wave combination condition, and output the resonance amplification factor;

[0227] The assessment report generation module is used to generate an engineering structure safety assessment report based on pile foundation stress parameters and resonance amplification factors.

[0228] In the above technical solution, the data acquisition module integrates multiple collection functions to comprehensively collect key data for port projects, laying a solid data foundation for subsequent analysis. The elevation model generation module generates a seafloor digital elevation model through spatial registration, improving the availability and accuracy of topographic data. The data processing module ensures the time consistency and reliable quality of time series data, enhancing the credibility of subsequent calculations and simulations. The model construction module constructs an unstructured grid and nested spatial discretization model, combined with appropriate discrete parameters, to accurately simulate complex hydrodynamic environments. The parameter optimization module utilizes a variational assimilation algorithm to optimize hydrodynamic parameters for more realistic simulations. The wave-current combined field generation module generates a dataset that comprehensively reflects wave-current interactions. The working condition construction module constructs a variety of tidal and wave combination working conditions, providing realistic scenarios for engineering structural analysis. The pile foundation force calculation module and the resonance factor calculation module accurately calculate relevant parameters. The assessment report generation module uses these parameters to generate scientific safety assessment reports, assisting in decision-making and optimization of port projects and ensuring project quality and safety.

[0229] Specifically, the data acquisition module uses a multi-beam bathymetry system to measure the seabed topography. Its multi-beam can obtain large-area seabed information and generate high-precision underwater point cloud data. Coastline point cloud data is collected through airborne lidar to quickly obtain three-dimensional coastline information over a large area. Anemometers and wind vanes are installed on observation towers at different heights to record multi-layer wind speed and direction data at 10-minute intervals. An acoustic Doppler current meter (ADCP) is installed on the survey ship to measure along the predetermined survey line, and the acoustic Doppler effect is used to obtain velocity profile data at different depths. Wave buoys monitor wave direction, wave height, period and other parameters in real time and transmit data wirelessly.

[0230] Elevation model generation module: Using geographic information system (GIS) software, based on the shared geographic coordinates of seabed topography data and shoreline point cloud data, a registration algorithm is adopted to unify the two into the same spatial coordinate system. For example, the data is imported into ArcGIS software, and after coordinate matching and data fusion, a seabed digital elevation model is generated.

[0231] Data Processing: Using Python's pandas and numpy libraries, and using Coordinated Universal Time (UTC) as the unified time base, the multi-layer wind speed and direction data, current profile data, and wave direction spectrum data are aligned to the same temporal resolution through methods such as temporal interpolation. Kalman filtering algorithms are then used to correct anomalies in the dataset, producing a quality verification dataset.

[0232] Model construction module: Based on the seafloor digital elevation model, the professional mesh generation software GAMBIT was used. The terrain gradient field was first extracted, and the mesh density control area was divided according to the gradient threshold. The Delaunay triangulation algorithm was used to generate the initial unstructured mesh. A secondary mesh was then refined using the frontier advancing method in the core reinforcement area. A three-level nested hierarchy was constructed, and the mesh resolution at each level was set. Discrete parameters were also set according to the requirements for wave directional spectrum discretization.

[0233] Parameter Optimization Module: A two-stage variational assimilation framework is constructed using numerical simulation software such as COMSOL Multiphysics or FVCOM. The tidal assimilation phase optimizes the bottom friction coefficient and vertical eddy viscosity coefficient using velocity profile data as constraints. The wave assimilation phase optimizes the JONSWAP peak enhancement factor and breakup index using wave directional spectrum data as constraints. A terrain correction cost function is designed and solved using the Quasi-Newton-L-BFGS algorithm. Dynamic weights are applied to conflicting regions of data, and the optimized parameter combination and confidence intervals are output.

[0234] The wave-current combined field generation module utilizes numerical simulation software to solve the wave-current coupled equations in wave-dominated, tidal-dominated, and combined correction phases. A turbulent kinetic energy compensation model is employed in the wave-breaking zone, and the correction term calculation step size is adjusted based on the terrain gradient. Bidirectional coupling iterations are performed to generate a combined wave-current field dataset.

[0235] The working condition construction module utilizes the boundary condition setting function of the numerical simulation software to set wave incidence boundary conditions based on the measured wave direction spectrum. Wave refraction correction angles are added in nearshore shallow waters based on the seafloor topography gradient. The main tidal components are selected, and a tidal level time series is generated based on the harmonic constant. The confidence level and terrain attenuation factor are considered to construct a combined tidal and wave working condition.

[0236] Pile foundation force calculation module: Import the wave-current combined field data set into structural mechanics analysis software such as ANSYS or ABAQUS, input data at a time step of 0.1 seconds, dynamically adjust the resistance coefficient and inertia coefficient based on the angle between the flow direction and the pile foundation axis and the vortex shedding zone, calculate the shear force and bending moment time history curves of each section of the pile foundation, and extract key indicators.

[0237] Resonance Factor Calculation Module: Utilizing engineering analysis software such as COMSOL Multiphysics, a three-dimensional characteristic equation is established based on the seafloor digital elevation model and the harbor contours to solve for the first six resonant modes. The dominant frequency of the wave spectrum is compared with the resonant mode frequencies to determine the degree of match. The resonance amplification factor is then output, taking into account the influence of secondary reflected waves at the breakwater opening.

[0238] Assessment report generation module: Use Python's pandas library to build a dual-index assessment matrix, divide the safety level into four levels, automatically generate an optimization suggestion library, and finally generate an engineering structure safety assessment report.

[0239] Although the embodiments of the present invention have been disclosed above, they are not limited to the applications listed in the description and implementation methods. They can be fully applied to various fields suitable for the present invention. For those familiar with the art, additional modifications can be easily implemented. Therefore, without departing from the general concept defined by the claims and the scope of equivalents, the present invention is not limited to the specific details and illustrations shown and described herein.

Claims

1. A tidal wave joint testing method based on multi-source data in port engineering construction, characterized by: The following steps are involved: S1. Acquire seabed topography data and shoreline point cloud data, collect multi-layer wind speed and direction data, current velocity profile data, and wave direction spectrum data; S2, spatially registering the terrain data with the shoreline point cloud data to generate a seabed digital elevation model; S3, perform spatiotemporal alignment processing on the multi-layer wind speed and direction data, velocity profile data, and wave direction spectrum data collected in S1 to generate a quality verification data set with temporal and spatial consistency; S4, constructing an unstructured grid based on the seabed digital elevation model generated in S2, generating a nested spatial discrete model, and setting wave directional spectrum discrete parameters in the nested spatial discrete model; S5, using the variational assimilation algorithm to optimize the quality verification data set of S3 and output the optimized hydrodynamic parameter combination; S6, based on the wave directional spectrum discrete parameters of S4 and the hydrodynamic parameters optimized by S5, solve the wave-current coupling equation and calculate the radiation stress gradient correction term to generate the wave-current joint field data set; S7. Configure the wave incident boundary conditions according to the wave direction spectrum data, synthesize the tide level boundary conditions according to the tidal harmonic constant, and combine the tide level and the design wave height to generate the tide wave combination condition; S8, calculate the pile foundation stress parameters based on the wave-current combined field data set generated in S6 and the tide-wave combined working condition in S7; S9, based on the wave-current joint field data set generated by S6 and the tide-wave combined working condition of S7, uses the resonance characteristic equation to calculate the frequency matching degree and output the resonance amplification factor; S10. Generate an engineering structure safety assessment report based on pile foundation stress parameters and resonance amplification factors.

2. The tidal wave joint testing method based on multi-source data in port engineering construction according to claim 1 is characterized in that: The step S4 specifically includes: S4.

1. Extract the terrain gradient field based on the seafloor digital elevation model, divide the grid density control area according to the gradient threshold, and set the core density area, transition area, and background sparse area; S4.

2. Generate an initial unstructured grid using the Delaunay triangulation algorithm. For the core densification area specified in S4.1, implement local secondary densification using the superimposed frontier advancing method, ensuring that the grid side lengths meet the following requirements: core area ≤ 2m, transition area ≤ 5m, and background area ≤ 10m. S4.

3. Construct a three-level nested hierarchy based on the unstructured grid generated in S4.2: a first-level global grid covering the entire computational domain with a resolution of 10-20m; a second-level local grid covering the port waters with a resolution of 5-8m; and a third-level refined grid focusing on the area surrounding the pile foundation cluster with a resolution of 1-3m. S4.4, coupling the nesting boundaries of the three nesting levels in S4.3 by a bidirectional flux interpolation method, and using virtual units to transfer the gradient flux of hydrodynamic parameters, wherein the hydrodynamic parameters include flow velocity and wave height; S4.

5. Based on the hydrodynamic parameter combination optimized in S5, dynamic adaptive encryption is triggered in the velocity shear layer and wave breaking zone, and the grid topology is adjusted in real time according to the grid hierarchy structure in S4.

3.

3. The tidal wave joint testing method based on multi-source data in port engineering construction according to claim 2 is characterized in that: The step S4 specifically further includes: S4.

6. Extract the shoreline main axis direction θ based on the seabed digital elevation model from S2, using the range of θ ± 30° as the wave incident main direction interval. Based on the multi-layer wind speed and direction data, current velocity profile data, and wave direction spectrum data collected from S1, discretize the direction spectrum at intervals of 5°-10° within the wave incident main direction interval, and at intervals of 15°-20° in other areas. S4.

7. Based on the wave directional spectrum energy distribution collected in S1, three frequency bands are divided: the core frequency band frequency interval Δf1 = 0.02 Hz, the transition frequency band Δf2 = 0.05 Hz, and the edge frequency band Δf3 = 0.1 Hz; S4.

8. Combine the three-level nested hierarchy constructed in S4.3 and configure the discrete parameters hierarchically in the three-level nested grids: first-level global grid, directional spacing of 15°, frequency spacing Δf3; second-level local grid, directional spacing of 10°, frequency spacing Δf2; third-level refined grid, directional spacing of 5°, frequency spacing Δf1; S4.9, introduces a wavenumber-water depth coupling constraint. Based on the seafloor digital elevation model in S2, when the local water depth h and wavelength L of the seafloor digital elevation model satisfy h / L≤0.05, trigger the diffraction enhancement mode of the directional spectrum, and improve the directional resolution of the area to 2°; S4.

10. Based on the velocity field assimilated in S5, the multi-directional wave interference compensation algorithm is activated for the area where the velocity-wave angle is greater than 60°, and the reverse secondary wave component is superimposed on the original directional spectrum discretized in S4.

6.

4. The tidal wave joint testing method based on multi-source data in port engineering construction according to claim 3 is characterized in that: The step S5 specifically includes: S5.

1. Construct a two-stage variational assimilation framework: In the tidal assimilation stage, the bottom friction coefficient and vertical eddy viscosity coefficient are optimized using the velocity profile data in the quality verification dataset generated by S3 as constraints; in the wave assimilation stage, the JONSWAP spectral peak enhancement factor and breakup index are optimized using the wave directional spectrum data in the quality verification dataset generated by S3 as constraints; S5.

2. Set a terrain correction cost function, which is a weighted sum of three terms: the mean square error term between the observed flow velocity and the model-calculated flow velocity, a coupling constraint term for the seafloor topography gradient and bottom shear stress, and a parameter sparsity regularization term. The weight coefficient of the topography gradient term is 0.3 to 0.6 times the weight of the observation term. S5.

3. Use the quasi-Newton-L-BFGS algorithm to solve the cost function in S5.

2. When the gradient descent rate for three consecutive iterations is less than the preset threshold, the step size is automatically reduced to 1 / 4 of the original value; S5.

4. For areas where velocity and wave data conflict, dynamic weights are applied according to the following rules: in the time dimension, the weight decays exponentially with a half-life of 2-3 hours, centered on the conflict start time; in the spatial dimension, the weight decays to 25%-30% of the center area within a 50-60 meter radius centered on the pile foundation. S5.

5. Output the optimized parameter combination and confidence interval. When the coefficient of variation of tidal parameters is greater than 15% or the coefficient of variation of wave parameters is greater than 25%, manual review is triggered.

5. The tidal wave joint testing method based on multi-source data in port engineering construction according to claim 4 is characterized in that: The step S6 specifically includes: S6.

1. Solve the wave-current coupling equations in stages: Based on the discrete parameters of the wave directional spectrum from S4, in the wave-dominated phase, the gentle slope equation is used to calculate the wave field, outputting wave height, wave direction, and breaking energy distribution. Based on the hydrodynamic parameters optimized in S5, in the tidal-dominated phase, the three-dimensional shallow water equations are solved to obtain the velocity and water level fields. In the joint correction phase, the wave radiation stress gradient is injected as a source term into the tidal current equation, and the velocity field is fed back into the wave refraction calculation. S6.

2. Calculate the radiation stress gradient correction term: Use a turbulent kinetic energy compensation model in the wave breaking zone to convert 10%-15% of the dissipated energy from breaking into radiation stress increments. Dynamically adjust the correction term calculation step size based on the terrain gradient of the seafloor digital elevation model in S2. The step size should be ≤2m in steep slopes and ≤10m in flat areas. S6.

3. Perform bidirectional coupling iterations: After each wave-current coupling calculation in S6.1, compare the relative error of wave height and velocity deviation between two consecutive iterations. Terminate the calculation if the threshold is not exceeded for three consecutive iterations. Otherwise, update the radiation stress field and iterate again based on the radiation stress gradient correction term in S6.

2. S6.

4. Generate a wave-current joint field dataset, which includes: a spatial vector field, including flow velocity, flow direction, wave height, and wave direction; an energy scalar field, including radiation stress intensity and turbulent kinetic energy density; and a dynamic characteristic field, including the location of the breaking zone and the vortex core area.

6. The tidal wave joint testing method based on multi-source data in port engineering construction according to claim 5 is characterized in that: The step S8 of pile foundation force calculation specifically includes: S8.

1. Input the wave-current combined field data set generated in S6 into the improved Morison equation at a time step of 0.1 seconds. Dynamically adjust the resistance coefficient C according to the angle between the flow direction and the pile axis and the turbulent kinetic energy in the wave-current combined field data of S6. d and the inertia coefficient C m : When the angle between the flow direction and the pile axis is greater than 45°, the resistance coefficient C is dynamically adjusted. d Increase by 10% to 15%; in the vortex shedding area, the inertia coefficient C m Reduce by 5% to 8%, and the vortex shedding zone is turbulent kinetic energy>0.3J / m 3 ; S8.

2. Using the coefficients adjusted in S8.1, calculate the shear force and bending moment time history curves for each section of the pile foundation, and then extract the maximum value, standard deviation, and fatigue damage index.

7. The tidal wave joint testing method based on multi-source data in port engineering construction according to claim 5 is characterized in that: The step S9 specifically includes: S9.

1. Establish a three-dimensional characteristic equation based on the seabed digital elevation model generated in S2 and the harbor contour line, and solve the first six resonance modes with a frequency range of 0.01 to 0.5 Hz; S9.

2. Compare the dominant frequency of the wave spectrum in the wave-current joint field dataset generated in S6 with the resonant modal frequency calculated in S9.1 and calculate the degree of match: when the dominant frequency difference is ≤0.02 Hz, it is determined to be a strong resonance risk, and the amplification factor is set to 1.5 to 3.0; when the dominant frequency difference is >0.1 Hz, it is determined to be no resonance, and the amplification factor is set to 1.0; S9.

3. Determine the location of the breakwater opening based on the seafloor digital elevation model generated in S2. Add a 10% to 20% increment to the amplification factor obtained in S9.2 to account for the effects of secondary reflections in this area.

8. A tidal wave combined testing system based on the method according to any one of claims 1 to 7, characterized in that: include: A data acquisition module is used to obtain seabed topography data, shoreline point cloud data, multi-layer wind speed and direction data, flow velocity profile data, and wave direction spectrum data; The elevation model generation module is used to spatially register terrain data with shoreline point cloud data to generate a seabed digital elevation model; A data processing module is used to perform spatiotemporal alignment processing on the collected multi-layer wind speed and direction data, flow velocity profile data, and wave direction spectrum data to generate a quality verification data set with temporal and spatial consistency; A model construction module is used to construct an unstructured grid based on the seabed digital elevation model, generate a nested spatial discrete model, and set wave directional spectrum discrete parameters in the nested spatial discrete model; Parameter optimization module, which is used to optimize the S3 quality verification data set using the variational assimilation algorithm and output the optimized hydrodynamic parameter combination; The wave-current combined field generation module is used to solve the wave-current coupling equation and calculate the radiation stress gradient correction term based on the combination of wave directional spectrum discrete parameters and hydrodynamic parameters to generate a wave-current combined field data set; The working condition construction module is used to configure the wave incident boundary conditions according to the wave direction spectrum data, synthesize the tide level boundary conditions according to the tidal harmonic constant, and combine the tide level with the design wave height to generate the tidal wave combination working condition; The pile foundation stress calculation module is used to calculate the pile foundation stress parameters based on the generated wave-current joint field data set and the tidal wave combined working condition; The resonance factor calculation module is used to calculate the frequency matching degree using the resonance characteristic equation based on the generated wave-current joint field data set and the tide-wave combination condition, and output the resonance amplification factor; The assessment report generation module is used to generate an engineering structure safety assessment report based on pile foundation stress parameters and resonance amplification factors.

Images