CN113033920A

CN113033920A - Method and system for predicting marine drift trajectory

Info

Publication number: CN113033920A
Application number: CN202110450309.0A
Authority: CN
Inventors: 张振昌; 曾银东; 李雪丁; 郭民权; 张少涵; 马锦山; 陈日清
Original assignee: FUJIAN MARINE FORECASTS; Fujian Agriculture and Forestry University
Current assignee: FUJIAN MARINE FORECASTS; Fujian Agriculture and Forestry University
Priority date: 2020-11-16
Filing date: 2021-04-25
Publication date: 2021-06-25

Abstract

The invention provides a method and a system for predicting a maritime drift trajectory, which solve the problem that the movement trajectory of a maritime drift object cannot be quickly and accurately determined in the prior art by acquiring the initial position of the drift object, determining the maritime parameter of a maritime area where the drift object is located according to the initial position and calculating the drift position of the drift object at a specific moment based on the maritime parameter, effectively shorten the maritime search and rescue time and improve the maritime search and rescue efficiency.

Description

Method and system for predicting marine drift trajectory

Technical Field

The invention relates to the field of object motion trajectory tracking, in particular to a method and a system for predicting an offshore drift trajectory.

Background

With the development of the sea by human being on a larger scale, the maritime activities such as marine transportation, marine fishery and the like are more and more frequent, the marine accident is more and more concerned by people, and the frequency of the maritime rescue or marine object search work is more and more increased. Time is the most valuable for maritime search and rescue operations. The search and rescue direction and the search and rescue area are accurate, and the efficiency of maritime search and rescue can be effectively improved.

Factors influencing the direction and the speed of the marine drifter in the sea are complex and intricate, and the drift of the marine drifter is mainly influenced by sea surface wind force, ocean current, gravity, buoyancy and waves. The gravity and the buoyancy are balanced, and the influence on the horizontal motion of the object can be ignored, so that the wind, the ocean current and the ocean wave are mainly acted on the object in the horizontal direction, the acting force of the wind and the ocean current is generated only when relative motion is generated between drift objects and the wind and the ocean current, and the waves have excitation effect on the drift objects.

The conventional method for predicting the marine drift trajectory usually considers one or more marine factors qualitatively and lacks quantitative calculation. Meanwhile, due to the calculation of the single drift object at the single position, the probability statistics is not carried out after a large number of particle simulation, and the uncertainty is large. Therefore, a method for predicting a marine drift trajectory is needed, which can comprehensively consider various marine factors, perform simulation statistics through a large number of particles, and accurately calculate the drift direction and speed of a marine drift object, thereby effectively shortening the time of marine search and rescue work, improving the working efficiency of marine search and rescue, and saving consumed manpower and material resources.

Disclosure of Invention

In view of the above, the present invention provides a method and a system for predicting a marine drift trajectory, so as to solve the problem that in the prior art, a motion trajectory of a marine drift object cannot be determined quickly and accurately according to the influence of various factors of the sea, thereby delaying a time for maritime search and rescue.

In order to solve the above technical problems, the proposed solution is as follows:

a marine drift trajectory prediction method comprises the following steps:

acquiring an initial position of a drift object;

determining ocean parameters of an ocean area where the drift objects are located according to the initial position;

calculating the drift position of the drift object at a specific moment based on the ocean parameters;

wherein the marine parameters include at least one of: turbulent whirl, surface flow field drift, wind-induced drift, three-dimensional flow field drift, or sea waves.

Preferably, the determining the ocean parameters of the ocean area where the drift object is located according to the initial position specifically includes:

and determining the random movement distance generated by turbulent vortex motion of the ocean area where the floater is located, the drift velocity of the ocean surface flow field and the wind-guiding drift velocity according to the initial position of the floater.

Preferably, the random movement distance Δ α of the drift generated by the turbulent vortex action is calculated by the following formula:

wherein, the delta alpha is the random movement distance in the alpha direction, and the alpha represents the x-axis, the y-axis or the z-axis direction; r is [ -1, 1]Uniformly distributed random number between, K_αIs a disturbance in the alpha directionThe coefficient, Δ t, is the time step.

Preferably, when the random movement distance of the drift object generated under the action of turbulent vortex motion is calculated, the disturbance coefficient K is determined by adopting a mode of simulating drift tracks by multiple particles_αThe values specifically include:

acquiring a plurality of particles with different initial positions, wherein the initial positions of the drifter are used as the circle center, and the preset distance is used as the radius;

respectively dividing the disturbance coefficient K_αSetting different values, and respectively calculating the drifting positions of the drifter and each particle under the action of turbulent vortex at the specific moment;

determining a confidence interval for the particle for a confidence level value based on a preset confidence level value and a plurality of positions of the particle;

judging the guarantee rate of the drifter in the confidence interval based on the position of the drifter;

determining the optimal disturbance coefficient K based on the assurance rate_αThe value is obtained.

Preferably, the drift velocity of the sea surface flow field is acquired based on ground wave radar observation data and/or ROMS three-dimensional flow field forecast.

Preferably, a four-dimensional variational data assimilation method is adopted, the ROMS three-dimensional flow field observation data are processed and converted into functional minimization with a dynamic mode as constraint, and the ROMS three-dimensional flow field forecast with the minimum deviation between the ROMS three-dimensional flow field observation data and the observation data in a specified time window is obtained by adjusting control variables.

Preferably, the wind guide drift velocity V_yIs calculated by the following formula:

V_y＝k·W

wherein W is the wind speed 10 meters above the sea surface, and k is the wind drag coefficient.

Preferably, the wind drag coefficient k is set based on the size of the sea wave.

Preferably, the calculating the drift position of the drift object at a specific moment based on the ocean parameters specifically includes:

calculating the drift position of the drift using the following formula:

wherein the content of the first and second substances,

the initial position of the drift, t0 is the initial time,

and the object drift speed at the time t is the composition of component speeds generated in the action process of a plurality of ocean parameters, and delta alpha is the random movement distance caused by turbulent vortex motion.

and calculating the drift position of the drift by adopting a MapReduce model.

Preferably, the calculating the drift position of the drift object by using the MapReduce model specifically includes:

recording predicted multi-particle information in a distributed file system (HDFS);

dividing the multi-particle information into M data segments, wherein M is an integer greater than 1;

distributing the M data fragments to different Worker workers in the MapReduce;

each Worker calculates the drift position of each particle according to the ocean parameters;

and summarizing the drift position of each particle and outputting a predicted drift track.

and determining the random movement distance generated by turbulent vortex motion of the ocean area where the floater is located and the drift caused by the three-dimensional flow field according to the initial position of the floater.

An offshore drift trajectory prediction system comprising:

the initial position acquisition device is used for acquiring the initial position of the drift;

the ocean parameter acquisition device is used for determining the ocean parameters of the ocean area where the drift objects are located according to the initial position;

the big data computing platform is used for computing the drifting position of the drifting object at a specific moment based on the ocean parameters;

An offshore drift trajectory prediction device comprising a memory, a processor and a computer program stored in said memory and executable on said processor, characterized in that said processor realizes the method steps of offshore drift trajectory prediction when executing said computer program.

A computer-readable storage medium, in which a computer program is stored which, when being executed by a processor, carries out the method steps of marine drift trajectory prediction.

According to the technical scheme, the method for predicting the marine drift track, provided by the embodiment of the application, solves the problem that the motion track of the marine drift object cannot be quickly and accurately determined in the prior art by acquiring the initial position of the drift object, determining the marine parameters of the marine area where the drift object is located according to the initial position and calculating the drift position of the drift object at a specific moment based on the marine parameters, effectively shortens the marine search and rescue time and improves the marine search and rescue efficiency.

Drawings

In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings used in the description of the embodiments or the prior art will be briefly described below, it is obvious that the drawings in the following description are only embodiments of the present invention, and for those skilled in the art, other drawings can be obtained according to the provided drawings without creative efforts.

Fig. 1 is one of the flow charts of the marine drift trajectory prediction method of the present invention.

Fig. 2 is a second flowchart of the marine drift trajectory prediction method of the present invention.

Fig. 3 is a schematic diagram illustrating an analysis of the influence of the wind conductance coefficient on the drift trajectory prediction according to the present invention.

FIG. 4 is a graph showing the density distribution and 95% confidence interval of the particles of the present invention.

FIG. 5 is a schematic diagram of the ground wave radar deployment and detection zone of the present invention.

FIG. 6 is a schematic diagram of the four-dimensional variational assimilation of the present invention.

FIG. 7 is a schematic flow chart of a big data platform drift trajectory prediction mode according to the present invention.

Fig. 8 is a schematic structural diagram of the marine drift trajectory prediction system of the present invention.

Fig. 9 is a schematic structural diagram of the marine drift trajectory prediction device of the present invention.

Detailed Description

The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.

The method for predicting the marine drift trajectory is suitable for the field of tracking the motion trajectory of an object.

The drift of the marine drifter is mainly influenced by the action of sea surface wind force, ocean current, gravity, buoyancy and waves. The gravity and the buoyancy are balanced, and the influence on the horizontal motion of an object can be ignored, so that the horizontal direction is mainly influenced by wind, ocean current and sea waves, the acting force of the wind and the ocean current is generated only when relative motion is generated between drift objects and the wind and the current, the action of waves on the drift objects is not completely clear, but the waves have an excitation effect on the drift objects on the water. It is generally believed that the effect of waves is negligible at wavelengths larger than the drift size.

As shown in fig. 1 and 2, the method for predicting an offshore drift trajectory provided by the present invention specifically includes:

step 101, obtaining an initial position of a drift object.

Assume an initial position of

Drift to a new position after a certain time step Δ t

And 102, determining ocean parameters of the ocean area where the drift objects are located according to the initial position.

The marine parameters mainly include the following parameters: turbulent whirl, surface flow field drift, wind-induced drift, three-dimensional flow field drift, or sea waves.

In one embodiment of the invention, the random movement distance generated by turbulent whirling of the ocean region in which the floater is located, the drift velocity of the ocean surface flow field and the wind-guiding drift velocity are determined according to the initial position of the floater.

The wind-induced drift refers to drift of drifted objects generated by dragging or wind-driven pressure, and the wind-driven pressure refers to air pressure on the water surface of the drifted objects in a certain motion state. Its wind-guiding drift velocity V_yThe calculation formula of (2) is as follows:

V_y＝k·W

in the formula, W is the wind speed 10 meters above the sea surface, k is the wind dragging coefficient, and the factors such as the drainage state, the immersion proportion of the drifted objects and the like are different along with the difference of the drifted objects. When the drifter is a ship or a buoy, k is between 0.01 and 0.05. Generally, 0.01-0.03 is taken, for example, 0.01 is taken for people falling into water, 0.03 is taken for life rafts, and 0.02 is taken for unpowered ships.

Here, the influence of sea waves needs to be considered, a drift object is easily surged by the sea waves in a place with large sea waves, the wind dragging coefficient at the moment can be increased accordingly, the value of k in the formula is modified correspondingly, and the wind dragging coefficient k is set based on the size of the sea waves.

The random movement distance delta alpha of the drift object under the turbulent vortex action is calculated by the following formula:

wherein, the delta alpha is the random movement distance in the alpha direction, and the alpha represents the x-axis, the y-axis or the z-axis direction; r is [ -1, 1]Uniformly distributed random number between, K_αAnd the disturbance coefficient in the alpha direction is shown, and the delta t is a time step.

As shown in fig. 3, when wind-induced drift and random motion are not considered, the water particle trajectory calculated by the small grid ocean current surface layer result is approximately consistent with the onsite floater trajectory, the drift trajectory trend of the floater is consistent with the motion trend of the Taiwan strait in summer circulation, the floater drifts from southwest to northeast, and meanwhile, the influence of tidal current reciprocating motion is obviously reflected, and the trajectory roundly moves in the north.

In order to compare and analyze simulation results under the influence of different wind conductance coefficients, the wind conductance coefficients of 0, 0.01, 0.02 and 0.03 are respectively selected in the test. The mean errors of the simulation and observation were 7.08, 6.26, 9.01 and 13.409 km, respectively. The comparison shows that the simulation effect is best when the wind conductance coefficient is 0.01.

In order to take into account the random motion of drifters under the turbulent vortex action, a multi-particle approach is used to simulate the many possibilities of drift trajectories of drifters. The randomness of turbulent eddies can be properly described using the monte carlo method. The monte carlo method, also called random sampling technique or statistical experiment method, can obtain the probability of occurrence of a certain event or the average value of a certain random variable by some experimental method when the problem to be solved is the probability of occurrence of the certain event or the expected value of the certain random variable. The random number used for perturbation is a uniform random number, i.e. the generated random number is a set of sample values of random variables distributed uniformly in the (-1, 1) interval. 10000 particles move independently under the action of turbulent vortex motion, and can cover a wider rangePossible trajectories of drifters. The method comprises the following steps: at the initial moment, 10000 particles are arranged in a circle by taking the site drift position as the center of a circle and taking 10km as the radius, the particles perform drift motion under the action of tidal current, ocean current, wind-induced drift and turbulent vortex, the positions of the particles are recorded, and the particles are stored once per hour. The results of each time were analyzed primarily as follows: a confidence level of 95% is set, that is, a minimum interval in which 9500 particles fall is set, as shown in fig. 4, in the case of 10000 particles, the region where 9500 particles with higher density are located is taken as an example, and whether the site point at the time is in the confidence region is judged. If the probability that the site point falls within the corresponding confidence region, i.e. the guarantee rate, is greater than 70% within the statistical time period, K_αIs reasonable, and, at K_αReasonable condition, the larger the guarantee rate is, the higher K is_αThe more reasonable, otherwise K_αIt is not reasonable.

Table 1 shows the various K_αIn this case, the number of times the site point falls in the model point (abbreviated as "number of fits" in the table) is counted in units of times, and the duration of the entire trajectory is 69 hours, once per hour. It can be seen that at K_αSimilarly, the higher the confidence level, the larger the confidence interval/search area involved, the more time the site point falls into the model point, and the same confidence level, K_αThe greater the number of times the site falls within the model confidence interval at 250. In general, the percentage of time that a field point falls within the confidence interval of a model point can only be met above 70% if the confidence level is 95%. The more times the site falls within the model confidence interval, the higher the likelihood of searching for missing objects in the future. K_αThe more reasonable the value of (a). It can be seen that K_αIt is reasonable to say that 250 has a high success rate compared with other values.

TABLE 1 different K_αStatistical results of probabilities of time-of-site points falling into model points

The drift velocity of the sea surface flow field is obtained based on the observation data of the ground wave radar and the ROMS three-dimensional flow field forecast.

The detection distance of the medium-range high-frequency ground wave radar developed by the national 863 plan can reach 200 kilometers. The radar system mainly comprises a compact antenna array (super-directional antenna technology), a high-stability receiving channel, a high-power all-solid-state transmitter and the like, and has ocean surface current, wave height and wind field inversion capability corresponding to the super-directional antenna technology. At present, the Fujian demonstration area is built in the Longhai and the east mountain, and can effectively monitor the sea surface wind field, the wave field and the flow field in the coastal areas of the southern Fujian, the shallow Taiwan, the Penghu water channel and the southwest Taiwan. The single station detection range is 200km, and one group of wind, wave and stream data is generated every 10 minutes.

As shown in fig. 5, the dragon sea station (LoHI) location: 118.14E 24.27N, east mountain station (DoSA) location: 117.48E 23.66N, with a radar normal of 132 degrees. The area where the vector flow is more reliable is the overlapping area which is 50-100 km from the radar and has an intersection angle less than 60 degrees with each normal line by comparing with the on-site measured data and other data. Wherein east mountain (DoSA) position: 117.48E, 23.66N; dragon sea (LoHI) location: 118.14E, 24.27N. The radar normal is 132 degrees.

The model used by the invention is developed by a Regional Ocean Modeling System (ROMS), the ROMS is a free surface, hydrostatic and three-dimensional nonlinear oblique pressure original equation model, and the model can simulate the movement of different scales, such as the circulation simulation of the global scale, the water level and flow field change caused by meteorological factors or astronomical tides of the mesoscale, and the fluid movement of small-scale river channels and the like can be calculated; an S coordinate system is adopted in the vertical direction, the vertical layering changes along with the terrain, and an orthogonal curve coordinate system and an Arakawa C differential grid are adopted in the horizontal direction; in order to improve the calculation efficiency, a mode separation method is adopted, the three-dimensional flow problem with the free surface is divided into a surface wave propagation problem (external mode) and an internal wave propagation problem (internal mode), and the external mode still uses explicit difference; the calculation of the swirl viscosity coefficient and the turbulence diffusion coefficient is determined by a 2.5-order turbulence model (Mellorand Yamada, 1982), which is free from the interference of human factors to a certain extent, and the model is widely applied to numerical simulation of different scales of circulation from a region to an ocean.

The four-dimensional variational method converts observation data processing into functional minimization taking a dynamic mode as constraint, and aims to minimize the deviation between a mode forecast result obtained by the control variable in a specified time window and actual observation data by adjusting the control variable. The four-dimensional variational method can take observation data of different moments, different regions and different properties, including unconventional observation data of satellites, radars and the like which are difficult to apply by an optimal interpolation method, into consideration as a whole, so that an initial field which is consistent with a forecasting mode is obtained.

The invention adopts a four-dimensional variational data assimilation method to process and convert the ROMS three-dimensional flow field observation data into functional minimization taking a dynamic mode as constraint, and obtains the ROMS three-dimensional flow field forecast with minimum deviation between the ROMS three-dimensional flow field observation data and the observation data in a specified time window by adjusting control variables.

Statistical basis for data assimilation: both the observed data and the model results have errors and some form of approximation, so they cannot be completely trusted during data assimilation, some way of compromise needs to be made to reduce the average error between the analysis field and the true value, and thus the uncertainty of the data needs to be considered quantitatively. The uncertainty can be analyzed by calibration or hypothesis of the error probability, and by means of probability theory, the problem can be transformed into an optimization problem, and it is hoped that the uncertainty can be reduced or not enlarged.

(1) Mode state vector and observation vector

A series of elements in the pattern used to describe the ocean state constitute a pattern state vector x, and the ocean true state vector and the pattern state variables are distinct, and the description of the true state vector tends to be more complex than that of the pattern state vector. We use x_tThe vector represents the true state of the sea at a certain moment in time, thisThe vector cannot be accurately obtained at any moment, and actually means the projection of a real state vector in a mode space in the current chapter; x is the number of_bThe vector represents the background state vector, which is a prior estimate of the true state prior to data assimilation, also known as the background field; forecasting the final analysis result required by the system, namely obtaining a mode state vector which is closer to the real state and is recorded as x by the system after data assimilation_aReferred to as the analysis field. The background field and the analysis field may be univariate or all mode predictor variables, and x is a three-dimensional distribution of (u, v, T, S, zeta, ubar, vbar). The mode variables are arranged in a grid and variables to form a vector of length n, where n is the product of the number of grids and the number of variables. Unknown 'true' x discretized by a pattern grid_tAlso a vector of length n. In the solution process, the analysis field x is often combined_aConverting to solve background field correction quantity delta x, wherein the correction quantity is called analysis increment and satisfies the following conditions:

x_a＝x_b+δx。

the target of the analysis is the desired analysis field x_aAs close to x as possible_t. In the analysis process, a series of observed values are used, and an observed vector y is formed by the observed values_oAssuming that the number of observed values is p, the observation vector y_oIs also p. To use these observation vectors in the analysis process, it is necessary to be able to compare with the state vectors. It is easy to view y if each dimension of the pattern space can be directly observed_oConsider the special case of a state vector, but in practical applications, the observations are always smaller than the mode variables, are not uniformly distributed in space, and the observation points are not necessarily at the grid points of the mode. Therefore, we need to introduce an observation operator H, which is a matrix with dimensions p × n, to transform the mode space to the observation space and then compare with the observation data. The H-operator is a combination of an interpolation operator, typically from a mode space to an observation space, and a transformation between a mode variable to an observation variable. The H matrix converts the vector of the pattern space into a corresponding value in the observation space, its transpose or the companion matrix H^TThe vector in the observation space is converted into a vector in the mode space.

The deviation between the observation and the mode vector is given by the following equation:

d＝y_o-H(x)

when using the background field x_bSubstituting the above equation, d is called the update variable; by x_aCarry-in is called the analytical residue.

(2) Error model

The best estimate of the sea state is obtained from a priori sea information (background field) and observations, and in order to optimally combine this information, statistical information about the errors of this information is needed to model the errors between these vectors and their true values. The observed value and the mode calculation result both contain errors, and the difference from the real state is marked as epsilon. Given that multiple analytical tests can be repeated under strictly identical conditions, the epsilon obtained in each test is not the same. Some statistics of epsilon, such as mathematical expectations and variances, can be derived from this experiment. At any one time of analysis, the specific error epsilon cannot be known, but at least the statistical characteristics of the error epsilon can be known, and a generally popular scalar density function model is a Gaussian function.

For practical data assimilation analysis, the dimensions of the background field and the analysis field are typically of the order of 10⁷The observation field is typically of the order of 10⁵Or 10⁶Therefore, the "error variance" must be replaced with an "error covariance matrix".

The error covariance matrix may be derived from a vector error

The product with its transpose is obtained and the expected value is obtained by averaging a number of cases:

the upper horizontal line here represents the mathematical expectation (i.e., the same as E ()). For a general error vector with the size of n, the covariance matrix is n multiplied by n matrix, which is a pairSymmetric and positive, with diagonal elements being the variance of the vector error component

The off-diagonal elements are the covariance between each pair of vectors. Removing each element by the product of the standard deviation

A correlation matrix can be obtained:

if it is

Is a variance diagonal matrix, then P ═ D^1/2CD^1/2。

For background field and observation field error estimation, the following concept is introduced:

background error: epsilon_b＝x_b-x_tThe covariance is B ═ E (ε)_bε_b ^T)。

And (3) analyzing errors: epsilon_a＝x_a-x_tThe covariance is A ═ E(s)_as_a ^T)。

And (3) observation error: epsilon_o＝y_o-h (x) with covariance R ═ E (epsilon)_oε_o ^T)。

At the same time, we assume that the observation and background errors are uncorrelated: e (ε)_oε_b ^T)＝0。

From the realization technology, the four-dimensional variation method is an important expansion of three-dimensional variation in time dimension, and the processing methods in the aspects of observation data processing, background field error physics, optimization algorithm and the like can be completely the same. The two are mainly different in that the four-dimensional variational system considers observation data and mode fields at a plurality of moments, so that a positive mode, a tangential mode and a concomitant mode need to be introduced into an assimilation system.

(2) Objective function

The four-dimensional variational method aims to enable the mode prediction to well fit the observed data in the time window, and because the mode prediction is obtained by a mode initial field forward integration mode, an optimal initial field is actually searched, so that the mode prediction can best fit the observed data in the assimilation time interval. The expression of four-dimensional variational assimilation by mathematical formulas is intended to minimize the contribution to the initial field x (t)₀) The objective function of (2):

wherein x₀The initial state of the prediction mode is,

is the background field at the initial moment in time,

is the i th observation data, H is the observation operator, x_iIs related to the observed data

The mode calculation result at the same time is obtained by the forward integration of the initial state, and can be recorded as

Indicating a pattern from t₀Integration of time to t_iThe result of the time of day.

The objective function is simply the initial state x for the mode₀And the analysis value at each time in the time interval is obtained from the initial state through forward mode integration, as shown in fig. 6. The objective function can be divided into two parts, one is to measure the deviation between the background field and the analysis field at the initial moment, and is marked as J_b(ii) a Another term measures the deviation between the observation and the analysis field, denoted J_oThe sum of each observed increment in the time interval.

The four-dimensional variational assimilation needs to solve a positive mode, a tangent mode and an adjoint mode in an iteration mode, and the calculation amount for directly solving the four-dimensional variational problem is huge. The objective function is constructed from the analysis variables, which ideally have the same resolution as the pattern variables, which, for high-resolution patterns, would undoubtedly make the minimization process very computationally expensive. Research shows that the target analysis variable is not solved directly, but the analysis increment is solved, and then the analysis increment is used for correcting the analysis variable, so that the calculation overhead can be greatly reduced. The basic idea is to express the objective function as an increment form, solve the analysis increment at low resolution, and then convert the analysis increment at low resolution into an analysis field at high resolution, wherein the analysis increment can be regarded as a disturbance of the control variable. Since the computation of the analysis increment can be performed at a lower resolution, the computation overhead of the four-dimensional variational method can be greatly reduced.

(3) Incremental form of an objective function

The increment method adopts an iteration process, takes a background field as an initial value of iteration at the beginning,

analyzing the values for the nth iteration

Is determined from the last analysis value

Plus the analysis increment for this iteration:

thus, the objective function can be rewritten as

Integrating the forecast in the model

In that

The treatment is obtained by using Taylor expansion:

ignoring its second order terms, there are:

in the preceding description, observation operator H has been mentioned_iIs a linear assumption, so that the operator H will be observed_iCarrying out linearization treatment:

substituting the above linearization process into the objective function can obtain the following incremental form:

wherein

This introduces a tangential mode L of the nonlinear mode M, whose effect is to obtain t_iDisturbance amount of time

Pattern computation

Starting from the initial state of the pattern, the integration is performed step by step along with the time step of the pattern, so that:

wherein

Represents from t_i-1Time analysis value t_iOne product step of the time of day analysis values.

Tangent linear mode

Represented by an initial time t₀Is a disturbance

Integral forward one step to get t_iDisturbance of time of day

If from time t₀To time t_iIn between, multiple steps are calculated, then the perturbation is from t₀Integral to t_iThe tangential pattern at a time is given by the product of the tangential patterns at each step. Namely, the following relations are provided:

solving the four-dimensional variational assimilation requires calculating the gradient of the objective function, namely the partial derivative of the objective function to the initial disturbance:

wherein

An adjoint pattern representing a tangential pattern, which can be written as a multiplication of an inverse integral since both the tangential and adjoint patterns are linearProduct:

in each iteration of the four-dimensional variation minimization process, the gradient is calculated, and in the forward integration process, the observation increment is calculated for each observation time

And is multiplied by

Followed by a companion mode

These weighted increments are integrated back to the initial time. After the gradient of the objective function is obtained, the four-dimensional variational method modifies the initial field state appropriately according to the gradient, then a new forward integration and a new observed increment are calculated again, and the process is repeated again.

The Hessian matrix of the objective function is:

in the iterative process, in order to accelerate the convergence speed of the minimized iterative process, the convergence speed can be accelerated by a preprocessing mode. Preconditioners are introduced to reduce the preconditioner number of the Hessian matrix of the objective function. The optimal preconditioner is to make the Hessian matrix an identity matrix, but since the observed terms of the Hessian matrix are too complex, B is generally used in practice^-1To determine the preconditioner:

δxⁿ＝Uvⁿ

wherein the transformation matrix U is Cholesky decomposition of the covariance matrix B of the background field error, that is, B ═ UU is satisfied^T. Introduction of a new variable vⁿObtaining a control variable v as a control variable for four-dimensional variational assimilationⁿThe objective function of (2):

the corresponding gradient is:

the corresponding Hessian matrix is:

after introducing the pre-condition treatment, the control variable is changed, the condition number is reduced, the inverse of solving the covariance matrix of the background field is avoided, the objective function is more spherical, and each iteration step can be closer to the minimum value of the objective function.

(4) Four-dimensional variational assimilation calculation process

Assuming that the tangential mode and the adjoint mode have been implemented, the specific computation steps for the four-dimensional variational are as follows:

1. calculating a nonlinear mode by using the initial field of the mode in an assimilation time window;

2. if the external circulation condition is met, the four-dimensional variation calculation is finished;

3. computing update vectors

4. Will control variable vⁿConversion to t₀Incremental analysis of time of day

5. The tangential mode is changed from t₀Integration of time to t_pTime of day, the analysis increment at each observation time is obtained

6. Transforming each analysis increment from the pattern space to the observation space by using an observation operator, and calculating

7. The calculation of the observation deviation, i.e. the addition of the analysis increment to the update vector, is performed in the observation space

And weighted multiplication

Matrices for transforming observed deviations into pattern space, i.e. for calculating observed deviations in pattern space

8. Taking the obtained observation deviation in the mode space as the input of the adjoint mode, and reversely integrating the adjoint mode;

9. adding the result obtained by the adjoint mode and the control variable to obtain the objective function gradient related to the control variable

10. Calculating a new control variable vⁿAnd when the convergence condition is not met, the step 4 is carried out to continue the inner loop iteration, otherwise, the step 1 is carried out to continue the outer loop iteration.

The calculation process of four-dimensional variation assimilation is divided into an outer loop and an inner loop, wherein the outer loop consists of an updating vector and an optimization part; the inner loop mainly completes the iteration step of the optimization algorithm, and solves the objective function and the gradient thereof through the forward integration of the tangential mode and the backward integration of the adjoint mode.

And 103, calculating the drift position of the drift object at a specific moment based on the ocean parameters.

Calculating the drift position of the drift using the following formula:

wherein the content of the first and second substances,

the initial position of the drift, t0 is the initial time,

This formula is called lagrangian tracking. Obviously, the core of the Lagrange tracking method is to solve

While

Is the total synthesis of the component speeds generated by the action process of each environmental dynamic factor. Considered by the present invention as a predictor of drift

The device mainly comprises a sea surface flow field containing tide and ocean current and a wind guide drifting speed, vertical movement can be considered in random movement, in a practical example, people falling into water do not drift on the sea surface, enter the lower part of the sea water and float on the sea surface after a period of time, the flow field inside the sea water is inconsistent with the flow field on the sea surface, and even the surface layer flow and the lower layer flow are opposite.

This constitutes a further embodiment of the invention, where the random movement distance generated by turbulent eddies in the ocean area where the float is located and the drift caused by the three-dimensional flow field are determined from the initial position of the float. Vertical random motion is considered during random motion, only bottom layer ocean current information needs to be considered after the seawater is in a lower layer, and information of a wind field and waves of a sea surface does not need to be used.

The prediction model is established on the basis of an ocean big data platform, a Hadoop distributed architecture is adopted, a distributed file system (HDFS) is constructed on the Hadoop distributed architecture, and a MapReduce calculation engine on the HDFS is operated. The Hadoop Distributed File System (HDFS) is designed to fit distributed file systems running on general purpose hardware (comfort hardware). HDFS is a highly fault tolerant system suitable for deployment on inexpensive machines. HDFS provides high throughput data access and is well suited for application on large-scale data sets. HDFS relaxes a portion of the POSIX constraints to achieve the goal of streaming file system data.

The drift trajectory prediction algorithm is designed to adopt a MapReduce calculation mode, which is a distributed big data programming model and is used for distributed operation of a large-scale data set. MapReduce is a simple software framework based on which written applications can run on large clusters of thousands of commercial machines and process data sets at the level of T in parallel in a reliable fault tolerant manner. The method supports parallel computing management of a mass data set, supports synchronous management of state, distribution, execution, exception and message of tasks, supports task execution process management, supports general MapReduce model execution, supports computing resource exception (computing node drop or network exception), fault tolerance of the computing tasks, dynamic expandability of the computing resources, dynamic expansion of the computing resources, and supports common compression algorithms for input and output of the computing tasks, such as GZIP, ZIP, BIZP2, LZO, SNAPPY and the like.

The principle of the programming model of MapReduce is as follows: an input set of key/value pairs is utilized to generate an output set of key/value pairs. Users of the MapReduce library express this calculation with two functions: map and Reduce. The user-defined Map function accepts an input key/value pair value and then generates a set of intermediate key/value pairs. The MapReduce library passes all intermediate value values with the same intermediate key value together to the Reduce function. The user-defined Reduce function accepts a set of intermediate key values and an associated value. The Reduce function combines these value values to form a smaller set of value values.

As shown in fig. 7, the calculation flow of MapReduce is as follows:

distributing the M data fragments to different Worker workers in the MapReduce;

Specifically, Map calls are distributed for execution on multiple machines by automatically splitting the input data for the Map call into a set of M data fragments. The incoming data segments can be processed in parallel on different machines. The intermediate key values generated by the Map call are split into R distinct partitions (e.g., hash (key) mod R) using a partition function, and the Reduce call is also distributed for execution on multiple machines. The number of partitions (R) and the partition function are specified by the user.

The MapReduce drift trajectory algorithm automatically divides input data (information of a plurality of particles) into a set of M data segments, and Map calls are distributed to a plurality of machines for execution. The incoming data segments can be processed in parallel on different machines. The intermediate key values generated by the Map call are split into R distinct partitions (e.g., hash (key) mod R) using a partition function, and the Reduce call is also distributed for execution on multiple machines. The number of partitions (R) and the partition function are specified by the user.

In summary, in order to consider disturbance caused by random motion to a drift trajectory, improve the success rate of maritime search and rescue, realize a Lagrange's method and Monte Carlo method combined application model algorithm on a big data platform for the first time, utilize the ability of the big data platform to expand a distributed computation framework, and simultaneously perform prediction computation on the drift trajectory of massive particles, cover the probability that drift objects may appear in a larger range, greatly improve the computation efficiency of a multi-particle probability prediction algorithm, and greatly reduce the computation time of the big data drift trajectory algorithm under the condition of the same computation ability; under the emergency condition, the computing resources can be linearly expanded, new distributed computing nodes can be rapidly deployed, the computing time can be reduced by times under the support of enough computing power, and the maritime search and rescue efficiency is improved.

Based on the same concept of the marine drift trajectory prediction method provided in the foregoing of the present invention, the present invention further provides a marine drift trajectory prediction system, as shown in fig. 8, the apparatus includes: the system comprises an initial position acquisition device 100, a marine parameter acquisition device 200 and a big data computing platform 300.

The initial position acquiring device 100 is used for acquiring the initial position of the drift.

And the ocean parameter acquisition device 200 is used for determining the ocean parameters of the ocean area where the drift objects are located according to the initial position.

And the big data computing platform 300 is used for computing the drift position of the drift object at a specific moment based on the ocean parameters.

Based on the same concept of the method for predicting the marine drift trajectory provided in the foregoing of the present invention, the present invention further provides a device for predicting the marine drift trajectory, as shown in fig. 9, the device includes: a memory 101, a processor 102, and a computer program stored in the memory and executable on the processor 102. The processor 102, when executing the computer program, performs the method steps of marine drift trajectory prediction.

Finally, it should also be noted that, in this document, the terms "comprises," "comprising," or any other variation thereof, are intended to cover a non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements does not include only those elements but may include other elements not expressly listed or inherent to such process, method, article, or apparatus. Without further limitation, an element defined by the phrase "comprising an … …" does not exclude the presence of other identical elements in a process, method, article, or apparatus that comprises the element.

The embodiments in the present description are described in a progressive manner, each embodiment focuses on differences from other embodiments, and the same and similar parts among the embodiments are referred to each other. For the device-like embodiment, since it is basically similar to the method embodiment, the description is simple, and for the relevant points, reference may be made to the partial description of the method embodiment.

The present application is described with reference to flowchart illustrations and/or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the application. It will be understood that each flow and/or block of the flow diagrams and/or block diagrams, and combinations of flows and/or blocks in the flow diagrams and/or block diagrams, can be implemented by computer program instructions. These computer program instructions may be provided to a processor of a general purpose computer, special purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, create means for implementing the functions specified in the flowchart flow or flows and/or block diagram block or blocks.

These computer program instructions may also be stored in a computer-readable memory that can direct a computer or other programmable data processing apparatus to function in a particular manner, such that the instructions stored in the computer-readable memory produce an article of manufacture including instruction means which implement the function specified in the flowchart flow or flows and/or block diagram block or blocks.

These computer program instructions may also be loaded onto a computer or other programmable data processing apparatus to cause a series of operational steps to be performed on the computer or other programmable apparatus to produce a computer implemented process such that the instructions which execute on the computer or other programmable apparatus provide steps for implementing the functions specified in the flowchart flow or flows and/or block diagram block or blocks.

The previous description of the disclosed embodiments is provided to enable any person skilled in the art to make or use the present application. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the generic principles defined herein may be applied to other embodiments without departing from the spirit or scope of the application. Thus, the present application is not intended to be limited to the embodiments shown herein but is to be accorded the widest scope consistent with the principles and novel features disclosed herein.

The flowchart and block diagrams in the figures illustrate the architecture, functionality, and operation of possible implementations of methods, apparatus, and computer program products according to various embodiments of the present invention. In this regard, each block in the flowchart and block diagrams may represent a module, segment, or portion of code, which comprises one or more computer-executable instructions for implementing the logical function(s). It should also be noted that, in some alternative implementations, the functions noted in the block may occur out of the order noted in the figures. It will also be noted that each block of the block diagrams and flowchart illustrations, and combinations of blocks in the block diagrams and flowchart illustrations, can be implemented by special purpose hardware-based systems that perform the specified functions or acts, or combinations of special purpose hardware and computer instructions.

Claims

1. A method for predicting an offshore drift trajectory, the method comprising:

acquiring an initial position of a drift object;

2. The method according to claim 1, wherein the determining the ocean parameters of the ocean area where the drift object is located according to the initial position specifically comprises:

3. The method of claim 2, wherein the random movement distance Δ α of the drifter under turbulent whirl is calculated by the following equation:

4. The method according to claim 3, wherein the disturbance coefficient K is determined by simulating a drift trajectory with a plurality of particles when calculating the random movement distance of the drift object under the action of turbulent vortex_αThe values specifically include:

5. The method of claim 2, wherein the surface flow field drift velocity is obtained based on ground wave radar observation data and/or a ROMS three-dimensional flow field forecast.

6. The method of claim 5, wherein the method of four-dimensional variational data assimilation is adopted to process and convert the ROMS three-dimensional flow field observation data into functional minimization using a dynamic mode as a constraint, and the ROMS three-dimensional flow field forecast with the minimum deviation from the observation data in a specified time window is obtained by adjusting control variables.

7. The method of claim 2, wherein the wind conduction drift velocity V_yIs calculated by the following formula:

V_y＝k·W

8. The method of claim 7, wherein the wind drag coefficient k is set based on a size of the ocean wave.

9. The method according to claim 2, wherein the calculating the drift position of the drift object at a specific moment based on the ocean parameters specifically comprises:

calculating the drift position of the drift using the following formula:

wherein the content of the first and second substances,

is the initial position of the drift, t₀As an initial moment of time, the time of day,

10. The method according to claim 9, wherein the calculating the drift position of the drift object at a specific moment based on the ocean parameters specifically comprises:

and calculating the drift position of the drift by adopting a MapReduce model.

11. The method according to claim 10, wherein the calculating the drift position of the drift object by using a MapReduce model specifically comprises:

distributing the M data fragments to different Worker workers in the MapReduce;

12. The method according to claim 1, wherein the determining the ocean parameters of the ocean area where the drift object is located according to the initial position specifically comprises:

13. An offshore drift trajectory prediction device, the device comprising:

14. An offshore drift trajectory prediction device comprising a memory, a processor and a computer program stored in said memory and executable on said processor, characterized in that said processor, when executing said computer program, carries out the steps of the method according to any one of claims 1 to 12.

15. A computer-readable storage medium, in which a computer program is stored which, when being executed by a processor, carries out the steps of the method according to any one of claims 1 to 12.