Trajectory optimization method and device for accurately deploying marine sensors under water

11686874 · 2023-06-27

Assignee

Inventors

Cpc classification

International classification

Abstract

The disclosure discloses a trajectory optimization method and device for accurately deploying marine sensors under water. The method includes the following steps. 1. Randomly select N sets of initial control variables within a range. 2. Input all of N sets of x.sub.i to the sensor's underwater glide kinematics and dynamics models, and calculate the smallest distance between N actual deployment positions and target deployment positions. 3. Determine whether the number of iterative operations is less than the preset value, if yes, perform global random walk and local random walk on N sets of x.sub.i, obtain N sets of x.sub.i again, and return to step 2; otherwise, go to step 4. 4. Output the control variable x.sub.i corresponding to Δs(x).sub.nminmin and the corresponding trajectory as the optimal control variable and optimal trajectory. The disclosure can improve the accuracy of prediction on the underwater three-dimensional trajectory of the marine sensor.

Claims

1. A trajectory optimization method for accurately deploying marine sensors under water, comprising the following steps: step 1: according to a preset search range x.sub.min≤x≤x.sub.max, N sets of values are randomly selected for initial control variables x=(U.sub.0, V.sub.0, W.sub.0, ω.sub.x0, ω.sub.y0, ω.sub.z0, ϕ.sub.0, θ.sub.0, ψ.sub.0), denoted as x.sub.i, wherein, i=1, 2 . . . , N, representing an i-th set, x.sub.min, x.sub.max are a lower limit and an upper limit of the search range, U.sub.0, V.sub.0 and W.sub.0 are components of an initial linear velocity in x-, y-, z-directions in an local coordinate system, ωx.sub.0, ω.sub.y0 and ω.sub.z0 are components of an initial angular velocity ω.sub.0 in the x-, y-, z-directions in the local coordinate system, and (ϕ.sub.0, θ.sub.0, ψ.sub.0) are initial Euler angles of the sensor; step 2: in an n-th iterative operation, all N sets of x.sub.i are input into an underwater glide kinematics and dynamics model for an underwater trajectory of the sensor, and N actual deployment positions are calculated, then, a minimum deployment error Δs(x.sub.i).sub.nmin between a target deployment position and the actual deployment position of each set of x.sub.i in the current iterative operation is calculated, thereafter, minimum distances corresponding to the N sets of x.sub.i are compared with each other, a minimum value is taken as a global minimum deployment error Δs(x).sub.nminmin of this iterative operation, wherein: when n=1, Δs(x.sub.i).sub.lmin=Δs(x.sub.i).sub.1, Δs(x).sub.1minmin is denoted as follows, and turn to step 3: Δs(x).sub.1minmin=Min{Δs(x.sub.1).sub.1, Δs(x.sub.2).sub.1, . . . Δs(x.sub.N).sub.1}; when n>1, Δs(x.sub.i).sub.nmin=Min{Δs(x.sub.i).sub.n, Δs(x.sub.i).sub.(n-1)min}, and then turn to step 3, under the circumstances, Δs(x).sub.nminmin=Min{Δs(x.sub.1).sub.nmin, Δs(x.sub.2).sub.nmin, . . . Δs(x.sub.N).sub.nmin} wherein a subscript n represents current iteration operation steps, Δs(x.sub.i).sub.n represents an error distance between the actual deployment position and the target deployment position in the n-th step of iterative operation in the i-th set, a subscript min is a mark for a minimum value in the set, and a subscript minmin is a mark for minimum value between the sets, Δs(x.sub.i).sub.nmin represents the internal minimum deployment error in the n-th step of iterative operation in the i-th set, and Δs(x).sub.nminmin represents a global minimum deployment error of a whole domain in the n-th step of iterative operation; step 3: determine whether n is less than a preset maximum number of iteration steps N.sub.0, if yes, a global random walk is performed, N sets of x.sub.i are updated, and then a local random walk is performed on the above updated N sets of x.sub.i according to a preset probability, the N sets of x are obtained again, and return to step 2; otherwise, go to step 4; step 4: the control variable x.sub.i corresponding to Δs(x).sub.nminmin and a corresponding trajectory are output as an optimal control variable and an optimal trajectory, wherein the optimal control variable and the optimal trajectory are used to improve an accuracy of a three-dimensional underwater trajectory prediction for the marine sensor to be deployed underwater.

2. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 1, wherein, in step 2, the underwater glide kinematics and dynamics model for the underwater trajectory of the sensor is as follows: U . = 1 m [ ( G - B ) sin θ + F dx + F e x ] - W ω y + V ω Z + L G ( ω y 2 + ω z 2 ) V . = 1 m [ - ( G - B ) cos θ sin ϕ + F dy + F ly + F e y ] - U ω z - L G ω . z W . = 1 m [ - ( G - B ) cos θ cos ϕ + F dz + F l z + F e z ] + U ω y + L G ω . y ω . y = 1 J y ( M l y + M dy + M e y + M G B y + m L G W . - m L G U ω y ) ω . z = 1 J z ( M l z + M dz + M e z - m L G V . - m L G U ω z ) in the formula, {dot over (U)}, {dot over (V)} and {dot over (W)} are acceleration components of the sensor in the x-, y-, z-directions in the local coordinate system; {dot over (ω)}.sub.y and {dot over (ω)}.sub.z are angular acceleration components of the sensor in the y-, z-directions in the local coordinate system; m is a mass of the sensor, U, V, and W are velocity components of the sensor in the x-, y-, z-directions in the local coordinate system, ω.sub.x, ω.sub.y and ω.sub.z are angular velocities of the sensor in the x-, y-, z-directions in the local coordinate system, ϕ, θ, and ψ are roll, pitch and yaw angles of the sensor in the local coordinate system, J.sub.y and J.sub.z are moment of inertia of the sensor with respect to the y-axis and z-axis, respectively, B is a buoyancy force of the sensor, M.sub.GBy is a hydrostatic moment caused by a relative offset L.sub.G between the buoyancy center B and a gravity center G, F.sub.e and M.sub.e are a hydrodynamic force and a moment of the sensor, F.sub.l and M.sub.l are a lifting force and a torque received by the sensor, F.sub.d and M.sub.d are a drag force and a torque received by the sensor, the subscripts x, y, z indicate that the corresponding parameter is a component on the corresponding coordinate axis; a relationship between the angular velocity and the Euler angles is obtained through a conversion rule for a local coordinate system and a global coordinate system: [ ω x ω y ω z ] = [ - ψ . sin ( θ ) ψ . sin ( ϕ ) cos ( θ ) ψ . cos ( ϕ ) cos ( θ ) ] + [ 0 θ . cos ( ϕ ) - θ . sin ( ϕ ) ] + [ ϕ . 0 0 ] wherein {dot over (ψ)}, {dot over (θ)}, {dot over (ϕ)} are change rates of a yaw angle ψ, a pitch angle θ and a roll angle ϕ that varies with time.

3. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 2, wherein the hydrodynamic forces in the y and z directions are calculated by the following expressions:
F.sub.ey={dot over (V)}m.sub.22−{dot over (ω)}.sub.ym.sub.24−UVM.sub.22(x.sub.t)+.sub.zx.sub.tM.sub.22(x.sub.t)
F.sub.ez={dot over (W)}m.sub.33−{dot over (ω)}.sub.zm.sub.35−UVM.sub.33(x.sub.t)+.sub.yx.sub.tM.sub.33(x.sub.t) hydrodynamic torques in the y and z directions are calculated according to the following formula:
M.sub.ey=−{dot over (W)}m.sub.35−{dot over (ω)}.sub.ym.sub.55−UW[m.sub.33+x.sub.tmM.sub.33(x.sub.t)]+.sub.y[m.sub.35−(x.sub.t).sup.2M.sub.33(x.sub.t)]
M.sub.ez=−{dot over (V)}m.sub.26−{dot over (ω)}.sub.zm.sub.66−UV[m.sub.22+x.sub.tM.sub.22(x.sub.t)]+.sub.z[m.sub.26−(x.sub.t).sup.2M.sub.22(x.sub.t)] in the formula, m.sub.22, m.sub.24, m.sub.26, m.sub.33, m.sub.35, m.sub.55 and m.sub.66 are three-dimensional additional mass coefficients of the sensor, and M.sub.22(x.sub.t) and M.sub.33(x.sub.t) respectively represent two-dimensional additional mass coefficients of an effective trailing edge in the y and z directions, x.sub.t is a longitudinal position of an effective trailing edge; for a sensor with a tip, M.sub.33(x.sub.t)=0, M.sub.22(x.sub.t)=0, for a marine sensor with two blunt ends, M.sub.33(x.sub.t)≠0, M.sub.22(x.sub.t)≠0.

4. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 2, wherein the drag force and the torque are generated by a relative translational movement of the sensor in contact with water, and an axial drag force F.sub.dx is calculated as follows: F dx = { - C 1 π D v ρ 2 L U .Math. "\[LeftBracketingBar]" U .Math. "\[RightBracketingBar]" - 1 8 ρ π C dt D 2 U .Math. "\[LeftBracketingBar]" U .Math. "\[RightBracketingBar]" , Laminar flow - C 2 1 2 ρ π D L U 2 - 1 8 ρ π C dt D 2 U .Math. "\[LeftBracketingBar]" U .Math. "\[RightBracketingBar]" , Transition flow C 3 1 2 ρ π D L U 2 - 1 8 ρ π C dt D 2 U .Math. "\[LeftBracketingBar]" U .Math. "\[RightBracketingBar]" , Turbulence wherein D is a cross-sectional diameter of the sensor, v is a kinematic viscosity of the water, ρ is a density of the water, C.sub.dt is a tangential shape drag coefficient of the sensor, C.sub.1, C.sub.2, and C.sub.3 are shape drag coefficients of the sensor when there are laminar flow, transition flow and turbulence, L is a length of the sensor; lateral drag forces F.sub.dy and F.sub.dz in the y and z directions are as calculated follows: F dy = 1 2 ρ C dy D - 0 5 L 0 . 5 L - ( V + ω z x ) .Math. "\[LeftBracketingBar]" V + ω z x .Math. "\[RightBracketingBar]" dx F d z = 1 2 ρ C d z D - 0 . 5 L 0 . 5 L - ( W - ω y x ) .Math. "\[LeftBracketingBar]" W - ω y x .Math. "\[RightBracketingBar]" dx resistance moments M.sub.dy and M.sub.dz caused by the lateral drag force are calculated as follows: M dy = 1 2 ρ C d z D - 0 . 5 L 0 . 5 L x ( W - ω y x ) .Math. "\[LeftBracketingBar]" W - ω y x .Math. "\[RightBracketingBar]" dx M d z = - 1 2 ρ C dy D - 0 . 5 L 0 . 5 L x ( V + ω z x ) .Math. "\[LeftBracketingBar]" V + ω z x .Math. "\[RightBracketingBar]" dx wherein C.sub.dy and C.sub.dz are respectively the drag coefficients in the y and z directions.

5. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 2, wherein if the marine sensor itself does not rotate around the axis, the lifting force F.sub.l and the moment M.sub.l are the lateral forces that are mainly caused by the asymmetry of a wake vortex distribution during the flow of water laterally crossing the sensor; if the sensor itself rotates, the water flow circulates laterally and generates pressure difference between two sides of the sensor, which forms the lifting force; the lifting force that is generated when the sensor does not rotate in the axial direction is caused by the asymmetry of the wake vortex distribution, and the lifting forces are calculated as follows: F l y = 1 2 ρ D - 0 . 5 L 0 . 5 L ( W - ω y x ) C l ( x ) dx F l z = - 1 2 ρ D - 0 5 L 0 . 5 L ( V + ω z x ) C l ( x ) dx M l y = 1 2 ρ D - 0 5 L 0 . 5 L x ( V + ω z x ) C l ( x ) dx M l z = 1 2 ρ D - 0 . 5 L 0.5 L x ( W - ω y x ) C l ( x ) d x wherein the lifting force coefficient C.sub.l(x) is calculated as follows: C l ( x ) = { K C L a sin [ π St ( t - 8 ) ] , t > 8 0 , t 8 , wherein St=f.sub.vD/V.sub.a is the Strouhal number, f.sub.v is a vortex shedding frequency, V.sub.a is a cross-flow velocity, K=±1 reflects a random behavior of the lifting direction during gliding, C.sub.La is a lift constant, and the time coefficient t′ is calculated as follows: t = s ( x ) R - 0 . 3 5 1 wherein s(x) is a displacement in the y direction when the water flows around a cross section of the sensor, and R is a radius of the cross section of the sensor; the lifting force and moment during axial rotation are calculated as follows: F ly = 1 2 ρ D W - 0 . 5 L 0 . 5 L C l ( x ) ( W - ω y x ) dx F lz = - 1 2 ρ D V - 0 5 L 0.5 L C l ( x ) ( V + ω z x ) dx M l y = 1 2 ρ D V - 0 . 5 L 0.5 L C l ( x ) ( V + ω z x ) xdx M l z = 1 2 ρ DC l ( x ) W - 0.5 L 0 . 5 L ( W - ω y x ) x d x wherein the lifting force coefficient C′.sub.l(x) is calculated as follows: C l ( x ) = { ω x D V a , ω x D V a 8 8 + 0.12 ( ω x D v a - 8 ) , ω x D V a > 8 .

6. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 1, wherein the sensor in the local body coordinate system is an axisymmetric slender body, and the change of an initial roll angle does not affect the movement of the slender body, therefore, the initial roll angle and an initial immersion depth are as follows:
ϕ.sub.0=0
d=−L/2 sin(θ.sub.0) thus, the initial control variable is simplified as x=(U.sub.0, V.sub.0, W.sub.0, ωx.sub.0, ω.sub.y0, ω.sub.z0, θ.sub.0, ψ.sub.0), and L represents the length of the sensor.

7. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 1, wherein in step 3, a global random walk is performed through Lévy flight, and global update is performed for N sets of control parameters x.sub.i, as shown in the following formula:
x.sub.i.sup.n+1=x.sub.i.sup.n+α.Math.Lévy(s,λ) in the formula, α=α.sub.0(x.sub.j.sup.t−x.sub.i.sup.t), α.sub.0 is a preset constant, x.sub.j represents the j-th set of control parameters x that is randomly generated, s represents the step length of the iterative operation, .Math. stands for item-based multiplication, and λ is an exponent of a Lévy distribution.

8. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 1, wherein after the global update in step 3, each set of control parameters x.sub.i has a probability of p.sub.a∈(0,1) to be locally update, and the update method of the local random walk is shown in the following formula:
x.sub.i.sup.n+1==x.sub.i.sup.n+1+ϵ.Math.(x.sub.p.sup.n+1−x.sub.q.sup.n+1), p,q=1,2, . . . N wherein x.sub.p and x.sub.q are two random solutions, ϵ conforms to an uniform distribution of an interval [0,1], and the superscript n+1 is the number of iterations in the next step.

9. A computer-readable/writable storage medium storing a computer program, wherein when the computer program is executed, the computer program performs the following steps: step 1: randomly generating values for N sets of initial control variables such as velocity, accelerations and orientations within a certain search range; inputting all of N sets of xi to an underwater glide kinematics and dynamics model for an underwater trajectory of the sensor, and calculating a smallest distance between N actual deployment positions and target deployment positions; step 2: inputting the selected N sets of initial control variables into the underwater glide kinematics and dynamics model for the underwater trajectory of the sensor; step 3: calculating the smallest distance between N actual deployment positions and target deployment positions; step 4: determining whether a number of iterative operations is larger than a pre-set value, if yes, outputting the trajectory and the corresponding control variable as an optimal trajectory and an optimal control variable, otherwise, performing global random walk and local random walk on the N sets of control variable, obtaining N sets of control variable again, and returning to step 2; step 5: save the optimal trajectory and the optimal control variable such as initial velocity, acceleration and orientation output in step 4 into the computer-readable/writable storage medium, wherein the optimal control variable and the optimal trajectory are used to improve an accuracy of a three-dimensional underwater trajectory prediction for the marine sensor to be deployed underwater.

10. A trajectory optimization device for accurately deploying a marine sensor under water, comprising: a computer-readable/writable storage medium of claim 9 and a processor, wherein the computer-readable/writable storage medium is able to store the computer program and the data output by an execution of the computer program, the processor is used for calling and processing computer programs stored in the computer-readable storage medium, once the computer program is called, it starts to solve the underwater glide kinematics and dynamics model for the underwater trajectory of the sensor, output the underwater trajectory of the sensor and the deployment distance error, finally, the optimal trajectory and optimal control variable such as initial velocity, acceleration and orientation can be obtained by sequential global and local random updates for control variable such as the initial velocity, acceleration and orientation.

11. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 2, wherein the sensor in the local body coordinate system is an axisymmetric slender body, and the change of an initial roll angle does not affect the movement of the slender body, therefore, the initial roll angle and an initial immersion depth are as follows:
ϕ.sub.0=0
d=−L/2 sin(θ.sub.0) thus, the initial control variable is simplified as x=(U.sub.0, V.sub.0, W.sub.0, ω.sub.x0, ω.sub.y0, ω.sub.z0, θ.sub.0, ψ.sub.0), and L represents the length of the sensor.

12. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 2, wherein in step 3, a global random walk is performed through Lévy flight, and global update is performed for N sets of control parameters x.sub.i, as shown in the following formula:
x.sub.i.sup.n+1=x.sub.i.sup.n+α.Math.Lévy(s,λ) in the formula, α=α.sub.0(x.sub.j.sup.t−x.sub.i.sup.t), α.sub.0 is a preset constant, x.sub.j represents the j-th set of control parameters x that is randomly generated, s represents the step length of the iterative operation, .Math. stands for item-based multiplication, and λ is an exponent of a Levy distribution.

13. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 2, wherein after the global update in step 3, each set of control parameters x.sub.i has a probability of p.sub.a∈(0,1) to be locally update, and the update method of the local random walk is shown in the following formula:
x.sub.i.sup.n+1=x.sub.i.sup.n+1ϵ.Math.(x.sub.p.sup.n+1−x.sub.q.sup.n+1), p,q=1,2, . . . N wherein x.sub.p and x.sub.q are two random solutions, ϵ conforms to an uniform distribution of an interval [0,1], and the superscript n+1 is the number of iterations in the next step.

14. A computer-readable storage medium, wherein the computer-readable storage medium stores a computer program, when the computer program is executed by a processor, the trajectory optimization method as claimed in claim 2 is implemented.

15. A trajectory optimization device for accurately deploying a marine sensor under water, comprising the computer-readable storage medium and the processor as claimed in claim 14, and the processor is configured for calling and processing computer programs stored in the computer-readable storage medium.

16. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 3, wherein the sensor in the local body coordinate system is an axisymmetric slender body, and the change of an initial roll angle does not affect the movement of the slender body, therefore, the initial roll angle and an initial immersion depth are as follows:
ϕ.sub.0=0
d=−L/2 sin(θ.sub.0) thus, the initial control variable is simplified as x=(U.sub.0, V.sub.0, W.sub.0, ω.sub.x0, ω.sub.y0, ω.sub.z0, θ.sub.0, ψ.sub.0), and L represents the length of the sensor.

17. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 3, wherein in step 3, a global random walk is performed through Lévy flight, and global update is performed for N sets of control parameters x.sub.i, as shown in the following formula:
x.sub.i.sup.n+1=x.sub.i.sup.n+α.Math.Lévy(s,λ) in the formula, α=α.sub.0(x.sub.j.sup.t−x.sub.i.sup.t), α.sub.0 is a preset constant, x.sub.j represents the j-th set of control parameters x that is randomly generated, s represents the step length of the iterative operation, .Math. stands for item-based multiplication, and λ is an exponent of a Lévy distribution.

18. The trajectory optimization method for accurately deploying marine sensors underwater according to claim 3, wherein after the global update in step 3, each set of control parameters x.sub.i has a probability of p.sub.a∈(0,1) to be locally update, and the update method of the local random walk is shown in the following formula:
x.sub.i.sup.n+1=x.sub.i.sup.n+1+ϵ.Math.(x.sub.p.sup.n+1−x.sub.q.sup.n+1), p,q=1,2, . . . N wherein x.sub.p and x.sub.q are two random solutions, ϵ conforms to an uniform distribution of an interval [0,1], and the superscript n+1 is the number of iterations in the next step.

19. A computer-readable storage medium, wherein the computer-readable storage medium stores a computer program, when the computer program is executed by a processor, the trajectory optimization method as claimed in claim 3 is implemented.

20. A trajectory optimization device for accurately deploying a marine sensor under water, comprising the computer-readable storage medium and the processor as claimed in claim 19, and the processor is configured for calling and processing computer programs stored in the computer-readable storage medium.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) FIG. 1 is the definition diagram of the coordinate system.

(2) FIG. 2 is a motion information diagram of the marine environment detection sensor at lateral dimension when not in axial rotation during the gliding process.

(3) FIG. 3 is a motion information diagram of the marine environment detection sensor at lateral dimension when in axial rotation during the gliding process.

(4) FIG. 4 is a display diagram of the underwater deployment operation of the marine environment detection sensor.

(5) FIG. 5 is a flow chart of the underwater trajectory optimization of the marine environment detection sensor.

(6) FIG. 6 is the motion trajectory diagram of the marine environment detection sensor after optimization under the target function of the minimum deployment error.

DESCRIPTION OF THE EMBODIMENTS

(7) In order to make the purpose, technical solutions, and advantages of the disclosure clearer, the following further describes the disclosure in detail with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described here are only used to explain the disclosure, but not to limit the disclosure. In addition, the technical features involved in the various embodiments of the disclosure described below can be combined with each other as long as they do not conflict with each other.

(8) As shown in FIG. 1 and FIG. 4, the marine environment detection sensor is deployed at a certain position underwater and performs a gliding motion until the sensor reaches a certain target position at the bottom of the water. In this embodiment, obtaining the minimum deployment error serves as an example, the underwater gliding motion trajectory of the marine environment detection sensor is optimized through the trajectory optimization method to set the main control variable value before the deployment. It should be noted that the trajectory optimization method is also applicable to the trajectory optimization of other marine sensors such as underwater seismic monitoring sensors, acoustic detection sensors, etc., and other underwater long-term observation type unpowered devices for underwater accurate deployment operations.

(9) First, the underwater glide kinematics model is established first. As shown in FIG. 1, the global coordinate system used as the inertial system is {I}-OXYZ, and the coordinate system is regarded as positive when the Z axis is at the top. In the meantime, the local coordinate system is represented by {B}-oxyz, which is a right-hand Cartesian coordinate system fixed at the gravity center of the marine environment detection sensor.

(10) In order to describe the motion of the sensor (in this embodiment, a slender cylinder) as a free fall in water, the following equations need to be solved in a local coordinate system.

(11) U . = 1 m [ ( G - B ) sin θ + F dx + F e x ] - W ω y + V ω z + L G ( ω y 2 + ω z 2 ) V . = 1 m [ - ( G - B ) cos θ sin ϕ + F dy + F l y + F e y ] - U ω z - L G ω . z W . = 1 m [ - ( G - B ) cos θ cos ϕ + F dz + F l z + F e z ] + U ω y + L G ω . y ω . y = 1 J y ( M l y + M dy + M e y + M G B y + m L G W . - m L G U ω y ) ω . z = 1 J z ( M l z + M dz + M e z - m L G V . - m L G U ω z )

(12) In the formula, {dot over (U)}, {dot over (V)} and {dot over (W)} are the acceleration components of the sensor in the x-, y-, z-directions in the local coordinate system; {dot over (ω)}.sub.y and {dot over (ω)}.sub.z are the angular acceleration components of the sensor in the y-, z-directions in the local coordinate system. m is the mass of the marine environment detection sensor, U, V, and W are the velocity components of the sensor in translation in the x-, y-, z-directions in the local coordinate system. ω.sub.x, ω.sub.y and ω.sub.z are the angular velocity of the sensor in the x-, y-, z-directions in the local coordinate system. ϕ, θ, ψ are the Euler angles of the sensor in the local coordinate system, which represent the axial rotation angle, pitch angle, and yaw angle, respectively. J.sub.y and J.sub.z are the moment of inertia of the sensor with respect to the y-axis and z-axis, respectively. Buoyancy B is the hydrostatic pressure item related to the drainage volume of the marine environment detection sensor. M.sub.GBy is the hydrostatic moment caused by the relative offset L.sub.G between the buoyancy center B and the gravity center G. F.sub.e and M.sub.e are the hydrodynamic force and moment of the sensor. F.sub.l and M.sub.l are the lifting force and torque received by the sensor. F.sub.d and M.sub.d are the drag force and torque received by the sensor. The subscripts x, y, z indicate that the corresponding parameter is the component on the corresponding coordinate axis.

(13) It should be noted that if the axial rotation rate along the x-axis is zero, that is, ω.sub.x=0 at the initial moment, the axial rotation rate is always very close to zero during the falling process. Otherwise, suppose ω.sub.x=c, and c is a fixed value throughout the process. In other words, in this embodiment, because the freely falling process is very short, for simplification, the possible deceleration of the axial rotation due to the viscosity of water is ignored.

(14) The relationship between the angular velocity of the marine environment detection sensor and the Euler angle is obtained through the conversion rule of the local coordinate system and the global coordinate system as follows:

(15) [ ω x ω y ω z ] = [ - ψ . sin ( θ ) ψ . sin ( ϕ ) cos ( θ ) ψ . cos ( ϕ ) cos ( θ ) ] + [ 0 θ . cos ( ϕ ) - θ . sin ( ϕ ) ] + [ ϕ . 0 0 ]

(16) A dynamic model is established for each force in the underwater gliding motion equation, and the kinematics equation and dynamics model are solved at the same time to obtain each transient motion parameter and magnitude of force. The underwater glide kinematics model and dynamics model can be established offline and then directly called online during the trajectory optimization process.

(17) (1) Calculation of Hydrodynamic Force and Moment:

(18) According to the potential flow theory under the slender body hypothesis, the hydrodynamic forces in the y and z directions are calculated as:
F.sub.ey=−{dot over (V)}m.sub.22−{dot over (ω)}.sub.ym.sub.24−UVM.sub.22(x.sub.t)+.sub.zx.sub.tM.sub.22(x.sub.t)
F.sub.ez={dot over (W)}m.sub.33−{dot over (ω)}.sub.zm.sub.35−UWM.sub.33(x.sub.t)+.sub.yx.sub.tM.sub.33(x.sub.t)

(19) The hydrodynamic torques in the y and z directions are calculated as follows:
M.sub.ey={dot over (W)}m.sub.35−{dot over (ω)}.sub.ym.sub.55−UW[m.sub.33+x.sub.tmM.sub.33(x.sub.t)]+.sub.y[m.sub.35−(x.sub.t).sup.2M.sub.33(x.sub.t)]
M.sub.ez=−{dot over (V)}m.sub.26−{dot over (ω)}.sub.zm.sub.66−UV[m.sub.22+x.sub.tM.sub.22(x.sub.t)]+.sub.z[m.sub.26−(x.sub.t).sup.2M.sub.22(x.sub.t)]

(20) In the formula, m.sub.22, m.sub.24, m.sub.26, m.sub.33, m.sub.35, m.sub.55 and m.sub.66 are the three-dimensional additional mass coefficients of the marine sensor, and M.sub.22(x.sub.t) and M.sub.33(x.sub.t) respectively represent the two-dimensional additional mass coefficients of the effective trailing edge in the y and z directions. x.sub.t is the longitudinal position of the effective trailing edge. For a marine environment detection sensor with a tip, the lateral force is zero: M.sub.33(x.sub.t)=0 M.sub.22(x.sub.t)=0. For a marine environment detection sensor with two blunt ends, M.sub.33(x.sub.t)≠0, M.sub.22(x.sub.t)≠0.

(21) (2) Calculation of axial drag force: since the water flow is generated axially, and is characterized in that its contact area with the marine environment detection sensor is large for a long time, and is mainly composed of frictional resistance and shape resistance, the axial drag force F.sub.dx is calculated as follows:

(22) F dx = { - C 1 π D v ρ 2 L U .Math. "\[LeftBracketingBar]" U .Math. "\[RightBracketingBar]" - 1 8 ρ π C dt D 2 U .Math. "\[LeftBracketingBar]" U .Math. "\[RightBracketingBar]" , Laminar flow - C 2 1 2 ρ π D L U 2 - 1 8 ρ π C dt D 2 U .Math. "\[LeftBracketingBar]" U .Math. "\[RightBracketingBar]" , Transition flow C 3 1 2 ρ π D L U 2 - 1 8 ρ π C dt D 2 U .Math. "\[LeftBracketingBar]" U .Math. "\[RightBracketingBar]" , Turbulence

(23) In the expression, D is the cross-sectional diameter of the marine environment detection sensor, v is the kinematic coefficient of water, ρ is the mass density of water, C.sub.dt is the tangential shape drag coefficient of the sensor, C.sub.1, C.sub.2, and C.sub.3 are the shape drag coefficient of the sensor when there are laminar flow, transition flow and turbulence, L is the length of the sensor.

(24) (3) Lateral drag force: Since water flow is generated in the manner of laterally surrounding the marine environment detection sensor, and is characterized in that, since the sensor has a large aspect ratio, it is estimated by the strip theory. The lateral drag force F.sub.dy and F.sub.dz are calculated as follows:

(25) F dy = 1 2 ρ C dy D - 0 . 5 L 0.5 L - ( V + ω z x ) .Math. "\[LeftBracketingBar]" V + ω z x .Math. "\[RightBracketingBar]" dx F dz = 1 2 ρ C dz D - 0 . 5 L 0.5 L - ( W - ω y x ) .Math. "\[LeftBracketingBar]" W - ω y x .Math. "\[RightBracketingBar]" dx

(26) Correspondingly, the resistance moments M.sub.dy and M.sub.dz that is caused by lateral drag force have the following expressions:

(27) M dy = 1 2 ρ C dz D - 0 . 5 L 0.5 L x ( W - ω y x ) .Math. "\[LeftBracketingBar]" W - ω y x .Math. "\[RightBracketingBar]" dx M dz = - 1 2 ρ C dy D - 0 . 5 L 0.5 L x ( V + ω z x ) .Math. "\[LeftBracketingBar]" V + ω z x .Math. "\[RightBracketingBar]" dx

(28) In the formula, C.sub.dy and C.sub.dz are the drag coefficients in the y and z directions, respectively.

(29) (4) If the marine environment detection sensor does not rotate axially, the lifting force and torque are the lateral forces that are generated mainly due to the wake vortex distribution asymmetry that is caused by the water flow that laterally surrounds the marine environment detection sensor. For the lifting force F.sub.ly and the lifting force F.sub.lz caused by the asymmetry of the wake vortex distribution, the lifting force moments M.sub.ly, and M.sub.lz are calculated as follows:

(30) F l y = 1 2 ρ D - 0 . 5 L 0.5 L ( W - ω y x ) C l ( x ) dx F l z = - 1 2 ρ D - 0 . 5 L 0.5 L ( V + ω z x ) C l ( x ) dx M l y = 1 2 ρ D - 0 . 5 L 0.5 L x ( V + ω z x ) C l ( x ) dx M l z = 1 2 ρ D - 0 . 5 L 0.5 L x ( W - ω y x ) C l ( x ) d x

(31) In the formula, the lifting force coefficient C.sub.1(x) is calculated as follows:

(32) C l ( x ) = { K C L a sin [ π St ( t - 8 ) ] , t > 8 0 , t 8

(33) In the expression, St=f.sub.vD/V.sub.a is the Strouhal number, f.sub.v is the vortex shedding frequency, V.sub.a is the cross-flow velocity. It should be noted that K=±1 reflects the random behavior of the lifting direction during gliding. In the deployment that is performed at different initial pitch angles as shown in FIG. 2, the lateral movement trajectories of the marine environment detection sensor are shown in the form of vibrations. C.sub.La is the lift constant, and the time coefficient t′ calculated is as follows.

(34) t = s ( x ) R - 0 . 3 5 1

(35) In the formula, s(x) is the displacement in the y direction when water flows around the cross section of the sensor, and R is the radius of the cross section of the sensor.

(36) (5) If the marine environment detection sensor rotates axially, the pressure difference between the two sides of the marine environment detection sensor is generated when the water flows laterally, which forms lifting forces F.sub.ly and F.sub.lz, and lift moments are M.sub.ly and M.sub.lz. FIG. 3 shows the changes in the three-dimensional motion trajectory of the marine environment detection sensor under different axial rotation rates.

(37) F l y = 1 2 ρ D W - 0 . 5 L 0 . 5 L C l ( x ) ( W - ω y x ) dx F l z = - 1 2 ρ D V - 0 . 5 L 0.5 L C l ( x ) ( V + ω z x ) dx M l y = 1 2 ρ D V - 0 . 5 L 0.5 L C l ( x ) ( V + ω z x ) xdx M l z = 1 2 ρ DC l ( x ) W - 0 . 5 L 0.5 L ( W - ω y x ) x d x

(38) In the formula, the lifting force coefficient C′.sub.l(x) is calculated as follows:

(39) 0 C l ( x ) = { ω x D V a , ω x D V a 8 8 + 0.12 ( ω x D v a - 8 ) , ω x D V a > 8

(40) FIG. 5 shows the process of optimizing the underwater trajectory of the marine environment detection sensor. First, based on the target function value output by the trajectory prediction model under the N sets of main initial control variables, combined with the global and local random walking algorithms, the main initial control variables are updated, and iterative operation is performed in sequence until the number of updates reaches the preset maximum number of steps. The specific steps are as follows:

(41) Step 1: The main initial control variables and initialization are determined. The initial control variables of the underwater three-dimensional trajectory of the marine environment detection sensor include 9 initial kinematics parameters x=(U.sub.0, V.sub.0, W.sub.0, ω.sub.x0, ω.sub.y0, ω.sub.z0, ϕ.sub.0, θ.sub.0, ψ.sub.0) and the initial position parameters (0, 0, d) of the marine environment detection sensor. U.sub.0, V.sub.0, W.sub.0 are the components of the initial linear velocity in the three coordinate axis directions in the local coordinate system, ω.sub.x0, ω.sub.y0, ω.sub.z0 are the components of the initial angular velocity in the three coordinate axis directions in the local coordinate system, ϕ.sub.0, θ.sub.0, ψ.sub.0 are the initial Euler angles, and d is the initial immersion depth.

(42) In this embodiment, since the slender body is axisymmetric, the change of the initial roll angle will not affect the movement of the slender body, so the initial roll angle and the initial immersion depth are as follows:
ϕ.sub.0=0
d=−L/2 sin(θ.sub.0)

(43) To sum up, in the scenario specified in this embodiment, first it is determined that the final eight main initial control variables of the marine environment detection sensor include:
x=(U.sub.0,V.sub.0,W.sub.0,ω.sub.x0,ω.sub.y0,ω.sub.z0,θ.sub.0,ψ.sub.0).

(44) Then, for specific task requirements, the initial control variables should satisfy a certain search range, as shown in the following formula:
x.sub.min≤x≤x.sub.max

(45) N sets of x.sub.i are selected randomly within the above range, wherein i=1, 2 . . . , N, then the initialization is completed in the first step.

(46) Step 2: All of N sets of x.sub.i are input into the underwater glide kinematics model and dynamics model of the sensor, and N actual deployment positions are calculated. Then the minimum distance between the target deployment position and the actual deployment position within each set are calculated, and the minimum distances in the N sets are compared with each other, and the minimum value is denoted as:
Δs(x).sub.nminmin=Min{Δs(x.sub.1).sub.nmin,Δs(x.sub.2).sub.nmin, . . . Δs(x.sub.N).sub.nmin}

(47) In the formula, the subscript n represents the current iteration operation steps, Δs (x.sub.i).sub.n represents the error distance between the actual deployment position and the target deployment position in the n-th step of iterative operation in the i-th set. The subscript min is the mark for minimum value in the set, and the subscript minmin is the mark for minimum value between the sets. Δs (x.sub.i).sub.nmin represents the internal minimum deployment error in the n-th step of iterative operation in the i-th set, and Δs(x).sub.nminmin represents the global minimum deployment error in the n-th step of iterative operation.

(48) When n=1, Δs(x.sub.i).sub.1 is directly used as Δs(x.sub.i).sub.1min, Δs(x.sub.i).sub.1minmin is solved, and then go to step 3;

(49) When n>1, Δs (x.sub.i).sub.nmin=Min{Δs (x.sub.i).sub.n, Δs (x.sub.i).sub.(n−1)min} is taken, then go to step 3.

(50) Step 2: The target function of the smallest deployment errors are obtained and compared. The global minimum deployment error is the minimum distance between the actual deployment position and the target deployment position, as shown in the following formula:
Δs(x).sub.nminmin=Min{Δs(x.sub.1).sub.nmin,Δs(x.sub.1).sub.nmin,Δs(x.sub.2).sub.nmin, . . . Δs(x.sub.N).sub.nmin}

(51) The N sets of x.sub.i are input into the underwater gliding motion and dynamics model to get the deployment error Δs (x.sub.i).sub.n calculated in the n-th step. For step n=1, the initialized N sets of x.sub.i are input into the underwater gliding motion and dynamics model to obtain the deployment error Δs(x.sub.i).sub.1 calculated in step 1, and the obtained Δs(x.sub.i).sub.1 is retained and used as the initial deployment error Δs(x.sub.i).sub.1min in the i-th set, N sets of Δs(x.sub.i).sub.1min are compared, and the minimum value Δs(x).sub.nminmin is taken as the initial deployment error.

(52) For step n>1, by comparing the Δs(x.sub.i).sub.n calculated in the n-step and the corresponding Δs(x.sub.i).sub.(n−1)min in the n−1-th step, the minimum value is selected as the internal minimum deployment error Δs(x.sub.i).sub.nmin in the i-th set in the n-step. N sets of Δs(x.sub.i).sub.nmin are compared, and the minimum value Δs(x).sub.nminmin is taken as the global minimum deployment error in the n-th step, which is denoted as:
Δs(x).sub.nminmin=Min{Δs(x.sub.1).sub.nmin,Δs(x.sub.2).sub.nmin, . . . Δs(x.sub.N).sub.nmin}
Specifically,
Δs(x.sub.i).sub.nmin=Min{Δs(x.sub.i).sub.n,Δs(x.sub.i).sub.(n−1)min}
Δs(x.sub.i).sub.n=√{square root over ((x.sub.t(x.sub.i).sub.n−a).sup.2+(Y.sub.t(x.sub.i).sub.n−b).sup.2+(Z.sub.t(x.sub.i).sub.n−c).sup.2)}

(53) (a,b,c) is the target deployment position, (x.sub.t(x.sub.i).sub.n, Y.sub.t(x.sub.i).sub.n, Z.sub.t(x.sub.i).sub.n) is the actual deployment position at the n-th step of iterative operation in the i-th set.

(54) Step 3: It is determined whether n is less than the preset maximum number of iteration steps N.sub.0. If n<N.sub.0, the N sets of initial control variables x in the n+1th step need to be updated through an optimization algorithm. The specific algorithm is as follows:

(55) (3.1) Global random walk algorithm. A global random walk is performed through Lévy flight, and the first update is performed for the N sets of control parameters x.sub.i, as shown in the following formula:
x.sub.i.sup.n+1=x.sub.i.sup.n+a.Math.Lévy(s,λ),

(56) In the formula, α=a.sub.0(x.sub.j.sup.n−x.sub.i.sup.n), and a.sub.0 is a constant. In most cases, α=O(1) can be adopted, that is, the order of a is usually 1. .Math. stands for item-based multiplication, and λ is the exponent of the Lévy distribution. The advantage of using the Lévy flight to perform a global random walk in this embodiment is that its step length is flexible, which can be large or small, free from the limitation of a single step length, and increases the search space per unit time. In other embodiments, other common global update algorithms can also be used instead.

(57) (3.2) Local random walk algorithm. After the first update (i.e., global update), there is a probability of p.sub.a ∈ (0, 1) for each set of control parameters x.sub.i in the N sets to be updated for the second time, mainly through a local random walk, as shown in the following formula:
x.sub.i.sup.n+1=x.sub.i.sup.n+1+ϵ.Math.(x.sub.p.sup.n+1−x.sub.q.sup.n+1), p,q=1,2, . . . N

(58) wherein p and q are two random solutions, and ϵ conforms to an uniform distribution in the interval [0,1].

(59) After completing the global update and partial update, return to step 2 and the iterative operation is performed in sequence until the number of updates n reaches the preset maximum number of times N.sub.0, that is, n≥N.sub.0, and then go to step 4.

(60) Step 4: The control variable x.sub.i corresponding to Δs (x.sub.i).sub.nmin, the internal minimum deployment error value Δs(x.sub.i).sub.nmin and the corresponding trajectory are output. FIG. 6 shows the optimized trajectory of the target function for the minimum deployment error in this embodiment.

(61) (6) A computer-readable/writable storage medium which stores a computer program. When the computer program is executed, the computer program performs the following steps:

(62) Step 1: randomly generate values for N sets of initial control variables such as velocity, accelerations and orientations within a certain search range. Input all of N sets of xi to the sensor's underwater glide kinematics and dynamics models, and calculate the smallest distance between N actual deployment positions and target deployment positions. 3. The disclosure can improve the accuracy of prediction on the underwater three-dimensional trajectory of the marine sensor.

(63) Step 2: input the selected N sets of initial control variables into the sensor's underwater glide kinematics model.

(64) Step 3: calculate the smallest distance between N actual deployment positions and target deployment positions.

(65) Step 4: determine whether the number of iterative operations is larger than the pre-set value, if yes, output the trajectory and the corresponding control variable as the optimal trajectory and optimal control variable. Otherwise, perform global random walk and local random walk on N sets of control variable, obtain N sets of control variable again, and return to step 2.

(66) Step 5: save the optimal trajectory and optimal control variable such as initial velocity, acceleration and orientation output in step 4 into the computer-readable/writable storage medium.

(67) (7) A trajectory optimization device for accurately deploying a marine sensor under water is provided, which includes the computer-readable/writable storage medium and a processor as described above, and the computer-readable/writable storage medium is able to store the computer program and the data output by the execution of the computer program, the processor is used for calling and processing computer programs stored in the computer-readable storage medium. Once the computer program is called, it starts to solve the sensor's underwater glide kinematics model, output the underwater trajectory of the sensor and the deployment distance error. Finally, the optimal trajectory and optimal control variable such as initial velocity, acceleration and orientation can be obtained by sequential global and local random updates for control variable such as the initial velocity, acceleration and orientation.

(68) In summary, the trajectory optimization method for accurate underwater deployment of marine environment detection sensors provided by the disclosure consists of a set of underwater three-dimensional trajectory prediction models and trajectory optimization algorithms. The underwater three-dimensional trajectory prediction model includes the underwater glide kinematics model and dynamics model of the marine environment detection sensor. The underwater glide kinematics model mainly includes the underwater gliding motion equation of the marine environment detection sensor and the conversion relationship between the local coordinate system and the linear velocity and angular velocity of the global coordinate system obtained through Euler angles. The dynamics model refers to the analysis and description of the force on the marine environment detection sensor when performing gliding motion under water. The trajectory optimization algorithm refers to the use of modern optimization algorithms to find the optimal input conditions and optimal trajectories that meet the underwater deployment operations. The input of the trajectory optimization algorithm is the main initial control variables of the three-dimensional trajectory prediction model, including the initial linear velocity, angular velocity and Euler angle at the moment when the marine environment detection sensor is released. The output of the algorithm is the trajectory of underwater gliding motion of the marine environment detection sensor, the deployment position error and so on.

(69) The underwater three-dimensional trajectory prediction model of the marine environment detection sensor and a trajectory optimization device provided by the disclosure can accurately describe and predict the underwater gliding motion of the marine environment detection sensor with six degrees of freedom. By taking into consideration the lateral force and displacement generated by the axial rotation or the asymmetric vortex, the accuracy of prediction on the underwater gliding trajectory of the marine environment detection sensor can be improved significantly. The trajectory optimization algorithm provided by the disclosure has flexible iteration steps, which can be large or small, free from the limitation of a single step length, so that the optimal solution search is efficient, and the trajectory can be optimized through multiple input variables.

(70) In addition, the method for optimizing the underwater trajectory of a marine environment detection sensor and a trajectory optimization device provided by the disclosure can be applied to trajectory optimization under the condition of a wider range of deployment operations, including high speed, short range, etc., and can be used to provide a more accurate guidance for optimal design and implementation of underwater deploying operation of marine environment detection sensors. In the meantime, the underwater trajectory optimization method can be used for performing accurate deployment and arrangement for other underwater long-term observation type unpowered devices such as underwater seismic monitoring, acoustic detection sensors, and so on.

(71) Those skilled in the art can easily understand that the above are only preferred embodiments of the disclosure and are not intended to limit the disclosure. Any modification, equivalent replacement and improvement, etc., made within the spirit and principle of the disclosure should be included in the protection scope of the disclosure.