CN111914411B - Full-attitude four-axis turntable frame angle instruction resolving method - Google Patents
Full-attitude four-axis turntable frame angle instruction resolving method Download PDFInfo
- Publication number
- CN111914411B CN111914411B CN202010697584.8A CN202010697584A CN111914411B CN 111914411 B CN111914411 B CN 111914411B CN 202010697584 A CN202010697584 A CN 202010697584A CN 111914411 B CN111914411 B CN 111914411B
- Authority
- CN
- China
- Prior art keywords
- attitude
- frame
- angle
- aircraft
- axis
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 78
- 238000004088 simulation Methods 0.000 claims abstract description 45
- 238000005457 optimization Methods 0.000 claims abstract description 35
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 29
- 239000011159 matrix material Substances 0.000 claims description 33
- 206010034719 Personality change Diseases 0.000 claims description 16
- 238000005096 rolling process Methods 0.000 claims description 13
- 230000008859 change Effects 0.000 claims description 12
- 238000004364 calculation method Methods 0.000 claims description 10
- 238000006243 chemical reaction Methods 0.000 claims description 9
- 238000013461 design Methods 0.000 claims description 6
- 238000005070 sampling Methods 0.000 claims description 5
- 230000005251 gamma ray Effects 0.000 claims description 4
- 230000008569 process Effects 0.000 claims description 3
- 230000004044 response Effects 0.000 claims description 3
- 238000010586 diagram Methods 0.000 description 9
- 230000036544 posture Effects 0.000 description 4
- 230000009466 transformation Effects 0.000 description 4
- 230000003044 adaptive effect Effects 0.000 description 2
- 238000004804 winding Methods 0.000 description 2
- 238000010521 absorption reaction Methods 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008878 coupling Effects 0.000 description 1
- 238000010168 coupling process Methods 0.000 description 1
- 238000005859 coupling reaction Methods 0.000 description 1
- 230000004069 differentiation Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000013178 mathematical model Methods 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/15—Vehicle, aircraft or watercraft design
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/04—Constraint-based CAD
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T10/00—Road transport of goods or passengers
- Y02T10/10—Internal combustion engine [ICE] based vehicles
- Y02T10/40—Engine management systems
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Geometry (AREA)
- General Physics & Mathematics (AREA)
- General Engineering & Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Aviation & Aerospace Engineering (AREA)
- Automation & Control Theory (AREA)
- Feedback Control In General (AREA)
- Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)
Abstract
The invention provides a method for calculating frame angle instructions of a full-attitude four-axis turntable, which comprises the steps of establishing a kinematics model of a frame angle of the four-axis turntable and a full-space aircraft attitude angle based on a double Euler method, initializing each frame angular position at a simulation initial moment, designing a turntable kinematics inverse solution algorithm based on a constraint optimization theory, calculating to obtain each frame angular position instruction of the four-axis turntable at the current moment according to the kinematics inverse solution algorithm, and designing a kinematics forward solution algorithm of the four-axis turntable to finally obtain all instructions of the frame angular positions. The method can ensure that the instruction resolving of the four-axis turntable does not generate large-amplitude resolving deviation jump in the full attitude on the premise of meeting the real-time performance, effectively avoiding the singular position of kinematics and the like, and eliminates the phenomenon of large-amplitude error caused by positive solution multi-resolution.
Description
Technical Field
The invention belongs to the technical field of semi-physical simulation of aircrafts, and particularly relates to a full-attitude four-axis turntable frame angle instruction resolving method.
Background
In aircraft semi-physical simulation and testing, a flight simulation turntable is one of the commonly used hardware devices. Among various simulation turntables, a three-axis simulation turntable is most commonly used, and can accurately simulate three attitude changes of a pitching attitude, a rolling attitude and a yawing attitude of an aircraft in space. However, due to the structural characteristics of the three-axis simulation turntable, the three-axis simulation turntable has a kinematic singular position and cannot adapt to the new requirement of full-attitude flight simulation. In order to avoid singular positions on the premise of ensuring the attitude of an aircraft in a simulation space, a four-axis turntable scheme is provided in engineering, and an additional secondary task is realized by increasing the redundant degree of freedom. However, unlike the three-axis turret, each frame angle corresponds to an attitude angle of the simulated aircraft one by one, and there is a coupling relationship between four frame angles and three attitude angles of the four-axis turret, that is, there may be countless sets of frame angles corresponding to a given set of attitude angles. For the calculation of each frame angle instruction of the four-axis turntable, two types of solutions mainly exist at present, namely a solution based on an adaptive adjustment algorithm and a solution based on an optimization problem containing constraint. In the first category, the solving scheme based on the adaptive adjustment algorithm solves the angular position instructions of the other three frames by actively designing the motion law of the base frame, so as to realize the simulation of the attitude of the aircraft. The scheme can effectively avoid the kinematic singular position of the four-axis turntable, but an ideal base frame motion law is difficult to design. In addition, when the full-attitude simulation of the aircraft is realized by the scheme, the situation that the angle position commands of each frame obtained by solving are discontinuous and the solving result has large-amplitude error jump can occur. The second category of solutions treats the instruction solving problem as a multi-solvability problem and converts it into an optimization problem to solve. According to the scheme, optimization variables, constraint conditions, objective functions, optimization algorithms and the like of the optimization problem are designed according to expected performance indexes of users, but all performance indexes such as resolving instantaneity, accuracy and avoidance of singular positions are difficult to consider simultaneously. In addition, when the scheme is adopted to realize the simulation of the full attitude of the aircraft, the phenomena of instruction discontinuity and large-amplitude error jump can also occur.
Disclosure of Invention
The invention aims to overcome the limitation that the existing scheme cannot simulate in a full-attitude range, and provides a full-attitude four-axis turntable frame angle instruction resolving method.
The invention is realized by the following technical scheme, and provides a full-attitude four-axis turntable frame angle instruction resolving method, which specifically comprises the following steps of:
the method comprises the following steps: according to the structure of a four-axis rotary table and the principle of simulating the change of the attitude of an aircraft, establishing a kinematics model of a frame angle of the four-axis rotary table and the attitude angle of the full-space aircraft based on a double Euler method, wherein the kinematics model comprises an angular position kinematics model and an angular velocity kinematics model;
step two: inputting aircraft attitude at simulation initial moment
Gamma (0) is the roll angle at the simulated initial moment,
measuring initial position phi of each frame angle of the four-axis turntable for simulating yaw angle at initial time and theta (0) for simulating pitch angle at initial time 0 ＝[φ 10 ,φ 20 ,φ 30 ,φ 40 ] T Initializing each frame angular position at the simulation initial moment according to the established angular position kinematic relationship between the four-axis turntable frame angle and the aircraft attitude angle to obtain phi (0) = [ phi ] (0) = 1 (0),φ 2 (0),φ 3 (0),φ 4 (0)] T ；
Step three: designing a turntable kinematics inverse solution algorithm based on a constraint optimization theory, and selecting an optimization variable in an optimization problem, namely each frame angular position phi = [ phi ] 1 ,φ 2 ,φ 3 ,φ 4 ] T Constraint conditions
m represents the number of constraint conditions, an objective function f (phi) and an optimization algorithm;
step four: input deviceAircraft attitude to be simulated at the current moment, namely the nth moment
Aircraft attitude simulation solution at the previous moment, namely the (n-1) th moment
And each frame angular position phi (n-1) = [ phi ] of the four-axis rotary table at the previous moment 1 (n-1),φ 2 (n-1),φ 3 (n-1),φ 4 (n-1)] T Calculating to obtain an angular position command phi (n) = [ phi ] of each frame of the four-axis rotary table at the current moment according to a kinematics inverse solution algorithm 1 (n),φ 2 (n),φ 3 (n),φ 4 (n)] T ；
Step five: designing a kinematics positive solution algorithm of the four-axis turntable, and according to the angular position command phi (n) = [ phi ] of each frame at the current moment 1 (n),φ 2 (n),φ 3 (n),φ 4 (n)] T And resolving to obtain the attitude simulation solution of the aircraft at the current moment
And will be
Inputting the motion data into inverse kinematics solution calculation at the next moment;
step six: and repeating the third step to the fifth step in an iterative manner, and resolving all the instructions of the angular position of the frame.
Further, the four-axis turntable is a vertical four-axis turntable.
Further, the first step specifically comprises:
a positive Euler method representation angular position kinematic model for describing the attitude of the aircraft by using 2-3-1 rotation sequence attitude angles is as follows:
R G b ＝R G l (1)
an inverse Euler method for describing the attitude of the aircraft by using the 2-1-3 rotation attitude angle represents an angular position kinematic model as follows:
R rG b ＝R G l (2)
wherein R is G b Representing an attitude matrix, R, describing the attitude change of the aircraft in order of 2-3-1 attitude angles rG b Representing an attitude matrix, R, describing the attitude changes of the aircraft in attitude angles of 2-1-3 rotation order G l Representing an attitude matrix for describing the attitude change of the vertical four-axis turntable load; gamma, the concentration of the gamma-rays,
theta respectively represents a rolling angle, a yaw angle and a pitch angle in the eulerian method; gamma ray r ，
θ r Respectively representing a rolling angle, a yaw angle and a pitch angle in an anti-Euler method; phi is a unit of 1 ，φ 2 ，φ 3 ，φ 4 Respectively showing the frame angles of four frames of the vertical four-axis turntable from outside to inside; g represents a ground coordinate system, b represents an aircraft coordinate system, and l represents a load coordinate system;
a positive Euler method representation angular velocity kinematic model for describing the attitude of the aircraft by using a 2-3-1 rotation sequence attitude angle is as follows:
an anti-Euler method for describing the attitude of the aircraft by using the 2-1-3 rotation attitude angle represents an angular velocity kinematics model as follows:
wherein, J G b The jacobian matrix is shown describing the aircraft as moving in 2-3-1 rotation order attitude angles,
respectively representing the roll angular velocity, the yaw angular velocity and the pitch angular velocity in a positive Euler method; j. the design is a square rG b The jacobian matrix is shown describing the aircraft as moving in 2-1-3 rotation order attitude angles,
respectively representing the rolling angular velocity, the yaw angular velocity and the pitch angular velocity in the anti-Eulerian method; j. the design is a square G l A jacobian matrix is shown describing the vertical four-axis turntable load motion,
respectively representing the frame angular velocities of the four frames of the vertical four-axis turntable from outside to inside; the vector of rotation of the aircraft is denoted by ω b ＝[ω bx ,ω by ,ω bz ] T Wherein ω is bx 、ω by 、ω bz Are respectively provided withThe projection of the rotation vector in the X-axis, Y-axis and Z-axis directions of an aircraft coordinate system; the vector of rotation of the load is ω l ＝[ω lx ,ω ly ,ω lz ] T Wherein ω is lx 、ω ly 、ω lz The projections of the rotation vector in the X-axis direction, the Y-axis direction and the Z-axis direction of a load coordinate system are respectively;
for each angular velocity in the formula (6) and the formula (7), angular position difference is adopted to approximate in engineering, the sampling time interval is set to be delta t, and the attitude angle of the aircraft at the current moment is respectively gamma n ，
θ n ，γ r n ，
θ r n The attitude angles at the previous time are gamma n-1 ，
θ n-1 ，γ r n-1 ，
θ r n-1 Then, there are:
let the attitude matrix R in the formulas (1) and (2) G b ＝R rG b = R, the (i, j) th entry in the matrix is denoted R ij Then, the conversion formula for calculating the inverse euler angle from the positive euler angle when the same attitude changes is described as follows:
when the positive Euler angle is solved by the inverse Euler angle, the conversion formula of the positive Euler angle is as follows:
in formulas (12) and (13), k =0 or 1;
defining a rounding function Δ (x) 1 ,x 2 ) As shown in equation (14):
Δ(x 1 ,x 2 )＝min(|x 1 -x 2 |,2π-|x 1 -x 2 |) (14)
when switching between positive and negative Euler angles in practical application, each Euler angle at the current moment is set as gamma n ，
θ n ，γ r n ，
θ r n The Euler angles at the previous time are gamma n-1 ，
θ n-1 ，γ r n-1 ，
θ r n-1 Introduction of a discriminant function d k
Respectively calculating d according to the values of k 0 And d 1 Get order d k And a group of attitude angles with small values are used as conversion results of positive and negative Euler angles at the current moment.
Further, in the second step, the first step,
aircraft attitude at initial time based on input simulation
And measuring the initial position phi of each frame angle of the four-axis rotary table 0 ＝[φ 10 ,φ 20 ,φ 30 ,φ 40 ] T Keeping the angular position phi of the base frame of the four-axis turntable at the initial simulation moment 1 (0)＝φ 10 At the same time, the other three frame angular positions phi 2 (0)，φ 3 (0)，φ 4 (0) Satisfy the requirement of
Wherein R is X (γ (0)) represents a rotation matrix rotating about the X-axis by an angle γ (0), R Z (theta (0)) represents a rotation matrix rotated by an angle of theta (0) about the Z axis,
indicating rotation about the Y axis
A rotation matrix of angles; r Y (φ 1 (0))，R X (φ 4 (0))，R Y (φ 3 (0))，R Z (φ 2 (0) For the same reason;
let R X (φ 4 (0))R Y (φ 3 (0))R Z (φ 2 (0) T) then there is
In the formula (17), S represents a sine function sin (·), and C represents a cosine function cos (·);
the initial angular positions of the pitching outer frame, the yawing middle frame and the rolling inner frame in the vertical four-axis turntable are provided
The solving formula is as follows:
wherein,
k is 1 or-1, two groups of actual attitude angles can be obtained through kinematics positive solution according to two groups of frame angular position instructions obtained through solution, and are respectively recorded as
And
the actual attitude angles of the aircraft at the current moment are gamma, phi and theta respectively; introducing a discriminant function d k
Respectively calculating d according to the values of k 1 ，d -1 Get order d k Group of small value
As the initial angular positions of the three frames inside the four-axis turntable at the current moment.
Further, in the third step,
the optimization variable of the optimization problem is selected as each frame angular position, and each frame angular position of the four-axis rotary table is respectively phi when a certain moment t is set 1 (t)，φ 2 (t)，φ 3 (t)，φ 4 (t); for the differential form of the frame angular positions involved in the optimization scheme, the differential form is adopted for replacement in the engineering, and the angular positions of the frames of the four-axis rotary table at the previous moment are respectively set asφ 1 (t-1)，φ 2 (t-1)，φ 3 (t-1)，φ 4 (t-1) the time interval between two moments is Δ t, the differential form of the angular position of the frames at t, i.e. the angular velocity of the frames
Can be finished to obtain:
the optimization problem objective function is designed as:
wherein, w 1 ，w 2 ，w 3 ，w 4 The weighting coefficients are respectively the angular velocity of each frame and represent the motion capability and the dynamic performance of each frame, and the four weighting coefficients meet w when being taken as values 1 ＞w 2 ＞w 3 ＞w 4 ＞0；w 1 ，w 2 ，w 4 The value of (a) is designed and set according to the dynamic response difference of the corresponding frame angle, and for the aircrafts with postures which are different and change within a large range, w 1 ＝w 2 ＝w 4 ＝Δt 2 (ii) a Weighting factor w of middle frame 3 The angular positions of the outer frame and the middle frame are jointly influenced, and the influences of the angular positions of the outer frame and the middle frame on the weighting coefficient of the middle frame are mutually independent; therefore, in engineering, in order to avoid the kinematic singular position of the four-axis turntable, the weighting coefficient w is used 3 Expressed as a two-dimensional normal distribution probability density function:
wherein x and y are respectively mod (| φ) 2 |,2π)，mod(|φ 3 L, 2 pi) representing the angular positions of the outer frame and the middle frame of the four-axis turntable; parameter K describes probability densityThe height of the degree function, parameter D describing the weighting coefficient of the intermediate frame when far from the singular position, μ 1 And mu 2 Describing the kinematic singular positions of the four-axis turntable, namely the corresponding angular positions of the outer frame and the middle frame; sigma 1 And sigma 2 Describing the change rate and the change amplitude of the frame weighting coefficient in the position close to the singular position;
the constraint conditions are equality constraint, and the equality constraint in engineering is acted by the full-attitude aircraft angular velocity kinematics equation based on the double Euler method, so that
Then, the equality constraint is:
wherein,
the optimization algorithm selects a Lagrange multiplier method, and Z represents an integer set.
Further, in step five, the designing of the kinematics forward solution algorithm of the four-axis turntable specifically includes:
according to four frame angular positions of the four-axis rotary table, the terminal load of the rotary table, namely three attitude angles of the equivalent simulated aircraft, satisfy the equation:
the solving formula of the load attitude angle is as follows:
in formula (29)
k takes 1 or-1, corresponding to two groups of results;
defining a rounding function Delta (x) 1 ,x 2 )＝min(|x 1 -x 2 |,2π-|x 1 -x 2 | in the kinematics positive solution process, each attitude angle of the load at a certain moment is set as gamma respectively n 、
θ n Then, each attitude angle at the previous time is γ n-1 、
θ n-1 Introducing a discriminant function d k
Wherein k is 1 or-1; respectively calculating d according to the values of k 1 ，d -1 Get an order d k And a set of attitude angles with small values is used as a result of the positive kinematic solution at the current moment.
The invention has the beneficial effects that:
1. the angular velocity kinematic relationship based on the double-Euler method established by the invention provides a theoretical basis of constraint conditions for the optimization problem converted by inverse kinematics solution in the third step, and successfully establishes the angular velocity kinematic relationship between the attitude angle of the full-attitude aircraft and each frame angle of the four-axis turntable. The established kinematic relationship between the attitude angle of the full-attitude aircraft and the angular velocity of each frame of the rotary table can ensure that a mathematical model which is full in attitude and can overcome the singularity problem of the Euler equation can be established for any given attitude angle form (such as 2-3-1 rotation sequence) with fixed rotation sequence and four-axis rotary tables (such as vertical or horizontal four-axis rotary tables) with different forms. In engineering, angular velocity is subjected to angular position differentiation to approximate substitution, and the real-time performance of calculation is guaranteed while the full-attitude simulation of the aircraft is eliminated.
2. The initialization algorithm of each frame angular position of the rotary table provided by the method successfully solves the problem that each frame angular position of the four-axis rotary table corresponding to the initial moment of the attitude simulation of the aircraft is not zero. During initialization, the base frame is adjusted with minimal amplitude to account for its worst dynamic performance. In addition, the initialization algorithm expands the value range of the inverse trigonometric function into the full space, and by designing the accept-and-reject function, the multi-solution problem caused by the expanded value range is rejected, so that the initialized instruction resolving can not generate large-amplitude resolving deviation jump in the full posture.
3. The kinematics positive solution improved algorithm provided by the method successfully solves the problems that the calculation result is greatly different from the theoretical value and the phenomenon of large-amplitude jump possibly caused by the traditional kinematics positive solution algorithm when the full-attitude simulation of the aircraft is carried out. The improved algorithm expands the value range of the inverse trigonometric function into the full space, and rejects the multi-solution problem caused by the expansion of the value range by designing a rejection function, thereby eliminating the phenomenon of large-amplitude error caused by the forward multi-solution.
4. The invention designs a method for dynamically changing weighting coefficients of a centering frame in a kinematic inverse solution scheme, and simultaneously considers the angular positions of an outer frame and a middle frame in a singular position configuration of a four-axis turntable. When the outer frame and the middle frame are close to the singular position, the weighting coefficient of the angular velocity of the middle frame is increased, the angular velocity of the middle frame is reduced, and the movement trend close to the singular position is restrained. When only one of the outer or middle frames approaches the singular position, the angular velocity of the middle frame does not change significantly. In engineering, the kinematic singular position of the four-axis turntable can be more simply and effectively avoided.
Drawings
FIG. 1 is a graph of a coordinate system transformation relationship based on 2-3-1 rotation Euler angles in the present invention;
FIG. 2 is a schematic structural diagram of a vertical four-axis turntable;
FIG. 3 is a schematic diagram of kinematic singular positions caused by structural limitations of a vertical four-axis turntable;
FIG. 4 is a graph of a coordinate system transformation relationship based on 2-1-3 rotation Euler angles in the present invention;
FIG. 5 is a schematic diagram illustrating respective range divisions of the dual-Euler method of the present invention;
FIG. 6 is a schematic diagram of a result of inverse kinematics solution in the present invention;
FIG. 7 is a flow chart of a method for calculating an angular position command of a frame of the full-attitude flight simulation turntable according to the invention;
FIG. 8 is a schematic view of an attitude and angular position command of an aircraft in an embodiment;
FIG. 9 is a schematic diagram of a result obtained by a four-axis turntable frame angular position calculation method in the embodiment;
FIG. 10 is a schematic diagram of a simulation solution of the attitude and angular positions of the aircraft obtained from the positive solution of the frame angles of the four-axis turntable in the embodiment;
FIG. 11 is a schematic diagram of a calculated deviation between an expected attitude angle of the aircraft and an attitude angle simulation solution in the embodiment.
Detailed Description
The technical solutions in the embodiments of the present invention will be described clearly and completely with reference to the accompanying 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.
With reference to fig. 1 to 11, the invention provides a method for calculating a frame angle instruction of a full-attitude four-axis turntable, which specifically comprises the following steps:
the method comprises the following steps: according to the structure of a four-axis turntable and the principle of simulating the attitude change of an aircraft, establishing a kinematics model of a frame angle of the four-axis turntable and an attitude angle of a full-space aircraft based on a double Euler method, wherein the kinematics model comprises an angular position kinematics model and an angular velocity kinematics model;
step two: inputting aircraft attitude at simulation initial moment
Gamma (0) is the roll angle at the simulated initial moment,
measuring initial position phi of each frame angle of the four-axis turntable to simulate yaw angle at initial moment, and theta (0) is pitch angle at the simulated initial moment 0 ＝[φ 10 ,φ 20 ,φ 30 ,φ 40 ] T Initializing each frame angular position at the simulation initial moment according to the established angular position kinematic relationship between the four-axis turntable frame angle and the aircraft attitude angle to obtain phi (0) = [ phi ] (0) = 1 (0),φ 2 (0),φ 3 (0),φ 4 (0)] T ；
Step three: designing a turntable kinematics inverse solution algorithm based on a constraint optimization theory, and selecting an optimization variable in an optimization problem, namely each frame angular position phi = [ phi ] 1 ,φ 2 ,φ 3 ,φ 4 ] T Constraint h (Φ) = [ h = [) 1 (Φ),…,h m (Φ)] T M represents the number of constraint conditions, an objective function f (phi) and an optimization algorithm;
step four: inputting the attitude of the aircraft to be simulated at the current moment, namely the nth moment
Aircraft attitude simulation solution at the previous moment, namely the (n-1) th moment
And each frame angular position phi (n-1) = [ phi ] of the four-axis rotary table at the previous moment 1 (n-1),φ 2 (n-1),φ 3 (n-1),φ 4 (n-1)] T Calculating to obtain an angular position command phi (n) = [ phi ] of each frame of the four-axis rotary table at the current moment according to a kinematics inverse solution algorithm 1 (n),φ 2 (n),φ 3 (n),φ 4 (n)] T ；
Step five: designing a kinematics positive solution algorithm of the four-axis turntable, and according to the angular position command phi (n) = [ phi ] of each frame at the current moment 1 (n),φ 2 (n),φ 3 (n),φ 4 (n)] T And resolving to obtain the attitude simulation solution of the aircraft at the current moment
And will be
Inputting the motion data into inverse kinematics solution calculation at the next moment;
step six: and repeating the third step to the fifth step in an iterative manner, and resolving all the instructions of the angular position of the frame.
The four-axis rotary table is a vertical four-axis rotary table. As shown in fig. 2, fig. 2 is a schematic structural diagram of a vertical four-axis turntable applicable to the four-axis turntable, in which yaw, pitch, yaw, and roll attitudes of an aircraft are simulated sequentially from the outside to the inside, and the aircraft rotates around the Y-axis direction, the Z-axis direction, the Y-axis direction, and the X-axis direction when the attitudes change. Wherein, OX B Y B Z B 、OX O Y O Z O 、OX M Y M Z M 、OX I Y I Z I Respectively representing a base frame coordinate system, an outer frame coordinate system, an intermediate frame coordinate system and an inner frame coordinate system. Fig. 3 shows the kinematic singular positions caused by the structural limitation of the vertical four-axis turntable. At this time, the angular positions of the outer frame and the middle frame of the turntable are phi 2 K pi, k is equal to Z and phi 3 K pi/2 + k pi, k epsilon Z. At the moment, the four-axis rotary table can only simulate the rolling and yawing or the rolling and pitching postures of the aircraft, and the pitching or yawing postures are lost, namely one degree of freedom is lost.
FIG. 1 is a graph of coordinate system transformation relation based on 2-3-1 rotation Euler angles in the present invention. Setting the ground coordinate system of the aircraft in the space as OX G Y G Z G The coordinate system of the projectile body is OX b Y b Z b The Euler angles in the sequence of 2-3-1 represent the attitude change of the aircraft from the position of the ground coordinate system to the position of the projectile coordinate system, and can be decomposed into successive OY windings G 、OZ’、OX 1 Three primitive rotations of the axis. FIG. 4 is a graph of coordinate system transformation relation based on 2-1-3 rotation Euler angles in the present invention. Set the ground coordinate system of the aircraft in the space as OX G Y G Z G The coordinate system of the projectile body is OX b Y b Z b The Euler angles in the sequence of 2-1-3 represent the attitude change of the aircraft from the position of the ground coordinate system to the position of the projectile coordinate system, and can be decomposed into successive OY windings G 、OX’、OZ 1 Three primitive rotations of the axis. The first step is specifically as follows:
the positive Euler method for describing the attitude of the aircraft by using the 2-3-1 rotation attitude angle represents an angular position kinematic model as follows:
R G b ＝R G l (1)
an inverse Euler method for describing the attitude of the aircraft by using the 2-1-3 rotation attitude angle represents an angular position kinematic model as follows:
wherein R is G b Representing an attitude matrix describing the attitude change of the aircraft in attitude angles of 2-3-1 turns, R rG b Representing an attitude matrix describing the attitude change of the aircraft in attitude angles of 2-1-3 turns, R G l Representing an attitude matrix for describing the attitude change of the vertical four-axis turntable load; the gamma-ray absorption material is a mixture of gamma,
theta respectively represents a rolling angle, a yaw angle and a pitch angle in the normal Eulerian method; gamma ray r ，
θ r Respectively representing a rolling angle, a yaw angle and a pitch angle in an anti-Euler method; phi is a unit of 1 ，φ 2 ，φ 3 ，φ 4 Respectively showing the frame angles of four frames of the vertical four-axis turntable from outside to inside; g represents a ground coordinate system, b represents an aircraft coordinate system, and l represents a load coordinate system;
a positive Euler method representation angular velocity kinematics model for describing the aircraft attitude by using a 2-3-1 rotation attitude angle is as follows:
an anti-Euler method for describing the attitude of the aircraft by using the 2-1-3 rotation attitude angle represents an angular velocity kinematics model as follows:
wherein, J G b The jacobian matrix is shown describing the aircraft as moving in attitude angles of 2-3-1 turns,
respectively representing the roll angular velocity, the yaw angular velocity and the pitch angular velocity in a positive Euler method; j is a unit of rG b The jacobian matrix is shown describing the aircraft as moving in 2-1-3 rotation order attitude angles,
respectively representing the rolling angular velocity, the yaw angular velocity and the pitch angular velocity in the anti-Eulerian method; j. the design is a square G l A jacobian matrix is shown describing the vertical four-axis turntable load motion,
respectively representing the frame angular velocities of four frames of the vertical four-axis turntable from outside to inside; the vector of rotation of the aircraft is denoted ω b ＝[ω bx ,ω by ,ω bz ] T Wherein ω is bx 、ω by 、ω bz Respectively projection of the rotation vector in X-axis, Y-axis and Z-axis directions of an aircraft coordinate system; the vector of rotation of the load is ω l ＝[ω lx ,ω ly ,ω lz ] T Wherein ω is lx 、ω ly 、ω lz Respectively at the load coordinate of the rotation vectorIs the projection in the X-axis, Y-axis and Z-axis directions;
for each angular velocity in the formula (6) and the formula (7), angular position difference is adopted to approximate in engineering, the sampling time interval is set to be delta t, and the attitude angle of the aircraft at the current moment is respectively gamma n ，
θ n ，γ r n ，
θ r n The attitude angles at the previous time are gamma n-1 ，
θ n-1 ，γ r n-1 ，
θ r n-1 Then, there are:
let the attitude matrix R in the formulas (1) and (2) G b ＝R rG b = R, item (i, j) in matrix is denoted R ij Then, the conversion formula for calculating the inverse euler angle from the positive euler angle when the same attitude changes is described as follows:
when the positive Euler angle is solved by the inverse Euler angle, the conversion formula of the positive Euler angle is as follows:
in formulas (12) and (13), k =0 or 1;
defining a rounding function Δ (x) 1 ,x 2 ) As shown in equation (14):
Δ(x 1 ,x 2 )＝min(|x 1 -x 2 |,2π-|x 1 -x 2 |) (14)
when switching between positive and negative Euler angles in practical application, each Euler angle at the current moment is set as gamma n ，
θ n ，γ r n ，
θ r n The Euler angles at the previous time are gamma n-1 ，
θ n-1 ，γ r n-1 ，
θ r n-1 Introducing a discriminant function d k
Respectively calculating d according to the value of k 0 And d 1 Get an order d k And a group of attitude angles with small values are used as conversion results of positive and negative Euler angles at the current moment.
FIG. 5 is a schematic diagram of the range division of the dual-Euler method according to the present invention. When the attitude angle of the 2-3-1 input is divided into two Euler angle regions, the second angle in the positive Euler angle, namely the pitching attitude angle position of the aircraft, is taken as a target when
Then, as shown by the blank part in fig. 5, the eueulerian method is applied. When the temperature is higher than the set temperature
Then, as shown by the shaded portion in FIG. 5, the inverse Euler method is applied.
In the second step, the first step is carried out,
according to the input aircraft attitude at the simulation initial moment
And measuring the initial position phi of each frame angle of the four-axis turntable 0 ＝[φ 10 ,φ 20 ,φ 30 ,φ 40 ] T Keeping the angular position phi of the base frame of the four-axis turntable at the initial simulation time 1 (0)＝φ 10 At the same time, the other three frame angular positions phi 2 (0)，φ 3 (0)，φ 4 (0) Satisfy the requirements of
Wherein R is X (γ (0)) represents a rotation matrix rotated about the X-axis by an angle of γ (0), R Z (theta (0)) represents a rotation matrix rotated by an angle theta (0) around the Z axis,
indicating rotation about the Y axis
A rotation matrix of angles; r is Y (φ 1 (0))，R X (φ 4 (0))，R Y (φ 3 (0))，R Z (φ 2 (0) The same applies;
let R be X (φ 4 (0))R Y (φ 3 (0))R Z (φ 2 (0) T) then there are
In the formula (17), S represents a sine function sin (·), and C represents a cosine function cos (·);
then a vertical four-shaft turntable is arranged inInitial angular positions of outer partial pitch frame, middle yaw frame, and inner roll frame
The solving formula is as follows:
wherein,
k is 1 or-1, corresponding to two groups of results. According to the two groups of frame angular position instructions obtained by resolving, two groups of actual attitude angles can be obtained through kinematics positive solution and are respectively recorded as
And
the actual attitude angles of the aircraft at the current moment are gamma, phi and theta respectively; introducing a discriminant function d k
Respectively calculating d according to the values of k 1 ，d -1 Get an order d k Group of small values
As the initial angular positions of the three frames inside the four-axis turntable at the current moment.
In the third step of the method, the first step,
in engineering, inverse kinematics solution problem is converted into optimization problem to solveAnd (6) determining. The optimization variable of the optimization problem is selected as each frame angular position, and when a certain moment t is set, each frame angular position of the four-axis rotary table is phi 1 (t)，φ 2 (t)，φ 3 (t)，φ 4 (t); for the differential form of the frame angular positions involved in the optimization scheme, the differential form is adopted for replacement in the engineering, and the angular positions of the frames of the four-axis rotary table at the previous moment are respectively set as phi 1 (t-1)，φ 2 (t-1)，φ 3 (t-1)，φ 4 (t-1), the time interval between the two moments being Δ t, the differential form of the angular position of the frames at time t, i.e. the angular speed of the frames
Can be finished to obtain:
the optimization problem objective function is designed as:
wherein, w 1 ，w 2 ，w 3 ，w 4 The weighting coefficients are the angular velocities of the frames respectively and represent the motion capability and the dynamic performance of the frames, and the four weighting coefficients meet w in value 1 ＞w 2 ＞w 3 ＞w 4 ＞0；w 1 ，w 2 ，w 4 The value of (a) is designed and set according to the dynamic response difference of the corresponding frame angle, and for the aircraft with the attitude which does not change at the same time and changes within a large range, generally, w is taken 1 ＝w 2 ＝w 4 ＝Δt 2 (ii) a Weighting coefficient w of middle frame 3 The angular positions of the outer frame and the middle frame are jointly influenced, and the influences of the angular positions of the outer frame and the middle frame on the weighting coefficient of the middle frame are mutually independent; therefore, in the engineering process, the air conditioner is provided with a fan,in order to avoid the singular position of the kinematics of the four-axis rotary table, a weighting coefficient w is used 3 Expressed as a two-dimensional normal distribution probability density function:
wherein x and y are respectively mod (| φ) 2 |,2π)，mod(|φ 3 L, 2 pi) representing the angular positions of the outer frame and the middle frame of the four-axis turntable; the parameter K describes the height of the probability density function, the parameter D describes the weighting coefficient of the intermediate frame, mu, away from the singular position 1 And mu 2 Describing the kinematic singular positions of the four-axis turntable, namely the corresponding angular positions of the outer frame and the middle frame; sigma 1 And sigma 2 Describing the change rate and the change amplitude of the frame weighting coefficient in the position close to the singular position; in engineering, each parameter is valued as mu 1 ＝π/2，μ 2 ＝π，K＝350，D＝1，σ 1 ＝π/9，σ 2 ＝π/9。
The constraint conditions are equality constraint, and the equality constraint in engineering is acted by the full-attitude aircraft angular velocity kinematics equation based on the double Euler method, so that
Then, the equality constraints are:
wherein,
the optimization algorithm selects a Lagrange multiplier method, and Z represents an integer set.
In the fifth step, the kinematics positive solution algorithm of the four-axis turntable is designed, and specifically comprises the following steps:
according to four frame angular positions of the four-axis rotary table, the terminal load of the rotary table, namely three attitude angles of the equivalent simulated aircraft, satisfy the equation:
the solving formula of the load attitude angle is as follows:
in formula (29)
k takes 1 or-1, corresponding to two groups of results;
defining a rounding function Delta (x) 1 ,x 2 )＝min(|x 1 -x 2 |,2π-|x 1 -x 2 | in the kinematics positive solution process, each attitude angle of the load at a certain moment is set as gamma respectively n 、
θ n Then, each attitude angle at the previous time is γ n-1 、
θ n-1 Introducing a discriminant function d k
Wherein k is 1 or-1; respectively calculating d according to the values of k 1 ，d -1 Get order d k And a group of attitude angles with small values are used as a kinematic positive solution result at the current moment.
Examples
In order to verify the effectiveness of the method, a vertical four-axis turntable is taken as a research object, the attitude of the aircraft to be simulated is input into a turntable system in a 2-3-1 rotation sequence, as shown in fig. 8, the sampling frequency is 2000Hz, namely the sampling time interval is delta t =0.0005s, and the total flight time of the aircraft is about 360s. The maximum pitch angle of the attitude input command of the aircraft reaches 90 degrees, and the condition of vertical launching of the aircraft can be simulated. The frame angle position instruction is resolved through a full-attitude four-axis rotary table frame angle instruction resolving method, and effectiveness of the method is verified. The method comprises the following specific steps:
(1) According to the flow shown in fig. 7, the angular positions of the frames of the four-axis turntable are initialized. All frames of the four-axis rotary table are in zero position before simulation begins, and the expected attitude angle of the aircraft is also in zero position at the initial moment of the simulation, so all frames of the four-axis rotary table obtained by initialization are still exactly kept in zero position, namely phi 1 (0)＝0，φ 2 (0)＝0，φ 3 (0)＝0，φ 4 (0)＝0，Φ(0)＝[0,0,0,0] T 。
(2) And inputting the expected aircraft attitude angle at each moment, the aircraft attitude simulation solution at the previous moment and the angular positions of the frames of the four-axis platform at the previous moment into the system, and resolving to obtain the angular position commands of the frames corresponding to the current moment.
(3) And (3) putting the angular position instructions of each frame at the current moment into kinematics positive solution calculation, and resolving to obtain a simulation solution of the attitude of the aircraft at the current moment.
(4) And (3) repeating the steps (2) and (3) according to the aircraft attitude command expected at the next moment, the aircraft attitude simulation solution obtained by calculation in the step (3) and each frame angular position command obtained by calculation in the step (2), and gradually obtaining all frame angular position commands.
According to the all-attitude four-axis turntable frame angle instruction resolving method, the final resolving result of each frame angle position is shown in fig. 9. The aircraft attitude simulation solution obtained through the kinematics positive solution is shown in fig. 10, and the resolving error between the attitude simulation solution and the aircraft expected attitude command is shown in fig. 11, so that the full-attitude four-axis turntable frame angle resolving method provided by the invention has good resolving accuracy, and meanwhile, large-amplitude error jump caused by full-attitude simulation can be avoided.
The method for calculating the frame angle instruction of the full-attitude four-axis turntable is described in detail, a specific example is applied to explain the principle and the implementation mode of the method, and the description of the embodiment is only used for helping to understand the method and the core idea of the method; meanwhile, for a person skilled in the art, according to the idea of the present invention, there may be variations in the specific embodiments and the application scope, and in summary, the content of the present specification should not be construed as a limitation to the present invention.
Claims (6)
1. A full-attitude four-axis turntable frame angle instruction resolving method is characterized by comprising the following steps: the method specifically comprises the following steps:
the method comprises the following steps: according to the structure of a four-axis rotary table and the principle of simulating the change of the attitude of an aircraft, establishing a kinematics model of a frame angle of the four-axis rotary table and the attitude angle of the full-space aircraft based on a double Euler method, wherein the kinematics model comprises an angular position kinematics model and an angular velocity kinematics model;
step two: inputting aircraft attitude at simulation initial moment
Gamma (0) is imitatedThe roll angle at the true initial time is,
measuring initial position phi of each frame angle of the four-axis turntable for simulating yaw angle at initial time and theta (0) for simulating pitch angle at initial time 0 ＝[φ 10 ,φ 20 ,φ 30 ,φ 40 ] T Initializing each frame angle position at the simulation initial moment according to the established angular position kinematic relationship between the four-axis turntable frame angle and the aircraft attitude angle to obtain phi (0) = [ phi ] (phi) 1 (0),φ 2 (0),φ 3 (0),φ 4 (0)] T ；
Step three: designing a turntable kinematics inverse solution algorithm based on a constraint optimization theory, and selecting an optimization variable in an optimization problem, namely each frame angular position phi = [ phi ] 1 ,φ 2 ,φ 3 ,φ 4 ] T Constraint h (Φ) = [ h = [) 1 (Φ),…,h m (Φ)] T M represents the number of constraint conditions, an objective function f (phi) and an optimization algorithm;
step four: inputting the attitude of the aircraft to be simulated at the current moment, namely the nth moment
Aircraft attitude simulation solution at the previous moment, namely the (n-1) th moment
And each frame angular position phi (n-1) = [ phi ] of the four-axis rotary table at the previous moment 1 (n-1),φ 2 (n-1),φ 3 (n-1),φ 4 (n-1)] T Calculating to obtain an angular position command phi (n) = [ phi ] of each frame of the four-axis rotary table at the current moment according to a kinematics inverse solution algorithm 1 (n),φ 2 (n),φ 3 (n),φ 4 (n)] T ；
Step five: designing a kinematics positive solution algorithm of the four-axis turntable, and according to the angular position command phi (n) = [ phi ] of each frame at the current moment 1 (n),φ 2 (n),φ 3 (n),φ 4 (n)] T And resolving to obtain the attitude simulation solution of the aircraft at the current moment
And will be
Inputting the motion data into inverse kinematics solution calculation at the next moment;
step six: and repeating the third step to the fifth step in an iterative manner, and resolving all the instructions of the angular position of the frame.
2. The method of claim 1, wherein: the four-axis rotary table is a vertical four-axis rotary table.
3. The method of claim 2, wherein: the first step is specifically as follows:
a positive Euler method representation angular position kinematic model for describing the attitude of the aircraft by using 2-3-1 rotation sequence attitude angles is as follows:
R G b ＝R G l (1)
an inverse Euler method representation angular position kinematic model for describing the attitude of the aircraft by using 2-1-3 rotation order attitude angles is as follows:
R rG b ＝R G l (2)
wherein R is G b Representing an attitude matrix, R, describing the attitude change of the aircraft in order of 2-3-1 attitude angles rG b Representing an attitude matrix, R, describing the attitude changes of the aircraft in attitude angles of 2-1-3 rotation order G l Representing an attitude matrix for describing the attitude change of the vertical four-axis turntable load; gamma, the concentration of the gamma-rays,
theta respectively represents a rolling angle, a yaw angle and a pitch angle in the normal Eulerian method; gamma ray r ，
θ r Respectively representing a rolling angle, a yaw angle and a pitch angle in the anti-Eulerian method; phi is a 1 ，φ 2 ，φ 3 ，φ 4 Respectively showing the frame angles of four frames of the vertical four-axis turntable from outside to inside; g represents a ground coordinate system, b represents an aircraft coordinate system, and l represents a load coordinate system;
a positive Euler method representation angular velocity kinematics model for describing the aircraft attitude by using a 2-3-1 rotation attitude angle is as follows:
an inverse Euler method representation angular velocity kinematics model for describing the attitude of the aircraft by using 2-1-3 rotation sequence attitude angles is as follows:
wherein, J G b The jacobian matrix is shown describing the aircraft as moving in 2-3-1 rotation order attitude angles,
respectively representing the roll angular velocity, the yaw angular velocity and the pitch angular velocity in a positive Euler method; j. the design is a square rG b The jacobian matrix is shown describing the aircraft as moving in 2-1-3 rotation order attitude angles,
respectively representing the roll angular velocity, the yaw angular velocity and the pitch angular velocity in an anti-Euler method; j is a unit of G l A jacobian matrix is shown describing the vertical four-axis turntable load motion,
respectively representing the frame angular velocities of the four frames of the vertical four-axis turntable from outside to inside; the vector of rotation of the aircraft is denoted by ω b ＝[ω bx ,ω by ,ω bz ] T Wherein ω is bx 、ω by 、ω bz Respectively projecting the rotation vector in X-axis, Y-axis and Z-axis directions of an aircraft coordinate system; the vector of rotation of the load is ω l ＝[ω lx ,ω ly ,ω lz ] T Wherein ω is lx 、ω ly 、ω lz The projections of the rotation vector in the X-axis direction, the Y-axis direction and the Z-axis direction of a load coordinate system are respectively;
for each angular velocity in the formula (6) and the formula (7), angular position difference is adopted to approximate in engineering, the sampling time interval is set to be delta t, and the attitude angle of the aircraft at the current moment is divided intoIs other than gamma n ，
θ n ，γ r n ，
θ r n The attitude angles at the previous time are gamma n-1 ，
θ n-1 ，γ r n-1 ，
θ r n-1 Then, there are:
let the attitude matrix R in the formulas (1) and (2) G b ＝R rG b = R, item (i, j) in matrix is denoted R ij Then, a conversion formula for solving the inverse euler angle from the positive euler angle when describing the same attitude change is as follows:
when the positive Euler angle is solved by the inverse Euler angle, the conversion formula of the positive Euler angle is as follows:
in formulas (12) and (13), k =0 or 1;
defining a rounding function Δ (x) 1 ,x 2 ) As shown in equation (14):
Δ(x 1 ,x 2 )＝min(|x 1 -x 2 |,2π-|x 1 -x 2 |) (14)
when switching between positive and negative Euler angles in practical application, each Euler angle at the current moment is set as gamma n ，
θ n ，γ r n ，
θ r n The Euler angles at the previous time are gamma n-1 ，
θ n-1 ，γ r n-1 ，
θ r n-1 Introducing a discriminant function d k
Respectively calculating d according to the value of k 0 And d 1 Get an order d k And a group of attitude angles with small values are used as conversion results of positive and negative Euler angles at the current moment.
4. The method of claim 3, wherein: in the second step of the method, the first step of the method,
aircraft attitude at initial time based on input simulation
And measuring the initial position phi of each frame angle of the four-axis rotary table 0 ＝[φ 10 ,φ 20 ,φ 30 ,φ 40 ] T Keeping the angular position phi of the base frame of the four-axis turntable at the initial simulation moment 1 (0)＝φ 10 At the same time, the other three frame angular positions phi 2 (0)，φ 3 (0)，φ 4 (0) Satisfy the requirement of
Wherein R is X (γ (0)) represents a rotation matrix rotating about the X-axis by an angle γ (0), R Z (theta (0)) represents a rotation matrix rotated by an angle theta (0) around the Z axis,
indicating rotation about the Y axis
A rotation matrix of angles; r Y (φ 1 (0))，R X (φ 4 (0))，R Y (φ 3 (0))，R Z (φ 2 (0) The same applies;
let R be X (φ 4 (0))R Y (φ 3 (0))R Z (φ 2 (0) T) then there is
In the formula (17), S represents a sine function sin (·), and C represents a cosine function cos (·);
the initial angular positions of the pitching outer frame, the yawing middle frame and the rolling inner frame in the vertical four-axis turntable are provided
The solving formula is as follows:
wherein,
k is 1 or-1, two groups of actual attitude angles can be obtained through kinematics positive solution according to two groups of frame angular position instructions obtained through solution, and are respectively recorded as
And
the actual attitude angles of the aircraft at the current moment are gamma, phi and theta respectively; introducing a discriminant function d k
Respectively calculating d according to the values of k 1 ，d -1 Get an order d k Group of small values
As the initial angular positions of the three frames inside the four-axis turntable at the current moment.
5. The method of claim 4, wherein: in the third step of the method, the first step,
the optimization variable of the optimization problem is selected as each frame angular position, and each frame angular position of the four-axis rotary table is respectively phi when a certain moment t is set 1 (t)，φ 2 (t)，φ 3 (t)，φ 4 (t); for the differential form of the angular position of the frame involved in the optimization scheme, the differential form is adopted in the engineering for replacement, and each frame of the four-axis rotary table at the previous moment is setAngular positions are respectively phi 1 (t-1)，φ 2 (t-1)，φ 3 (t-1)，φ 4 (t-1) the time interval between two moments is Δ t, the differential form of the angular position of the frames at t, i.e. the angular velocity of the frames
Can be finished to obtain:
the optimization problem objective function is designed as:
wherein w 1 ，w 2 ，w 3 ，w 4 The weighting coefficients are respectively the angular velocity of each frame and represent the motion capability and the dynamic performance of each frame, and the four weighting coefficients meet w when being taken as values 1 ＞w 2 ＞w 3 ＞w 4 ＞0；w 1 ，w 2 ，w 4 The value of (a) is designed and set according to the dynamic response difference of the corresponding frame angle, and w is the value of the aircraft with the attitude which is not changed at the same time and is changed in a large range 1 ＝w 2 ＝w 4 ＝Δt 2 (ii) a Weighting factor w of middle frame 3 The angular positions of the outer frame and the middle frame are influenced together, and the influence of the angular positions of the outer frame and the middle frame on the weighting coefficient of the middle frame is mutually independent; therefore, in the engineering, in order to avoid the singular position of the kinematics of the four-axis turntable, the weighting coefficient w is used 3 Expressed as a two-dimensional normal distribution probability density function:
wherein x and y are respectively mod (| φ) 2 |,2π)，mod(|φ 3 |,2π)，The angular positions of the outer frame and the middle frame of the four-axis turntable are shown; the parameter K describes the height of the probability density function, the parameter D describes the weighting coefficient of the intermediate frame, mu, away from the singular position 1 And mu 2 Describing the kinematic singular positions of the four-axis turntable, namely the corresponding angular positions of the outer frame and the middle frame; sigma 1 And sigma 2 Describing the change rate and the change amplitude of the frame weighting coefficient in the position close to the singular position;
the constraint conditions are equality constraint, and the equality constraint in engineering is acted by the full-attitude aircraft angular velocity kinematics equation based on the double Euler method, so that
Then, the equality constraint is:
wherein,
the optimization algorithm selects a Lagrange multiplier method, and Z represents an integer set.
6. The method of claim 5, wherein: in the fifth step, the kinematics positive solution algorithm for the four-axis turntable is designed, and specifically comprises the following steps:
according to four frame angular positions of the four-axis rotary table, the terminal load of the rotary table, namely three attitude angles of the equivalently simulated aircraft, satisfy the equation:
the solving formula of the load attitude angle is as follows:
in formula (29)
k takes 1 or-1, corresponding to two groups of results;
defining a rounding function Δ (x) 1 ,x 2 )＝min(|x 1 -x 2 |,2π-|x 1 -x 2 |) and in the process of kinematics positive solution, each attitude angle of the load at a certain moment is set as gamma respectively n 、
θ n Then, each attitude angle at the previous time is γ n-1 、
θ n-1 Introducing a discriminant function d k
Wherein k is 1 or-1; respectively calculating d according to the values of k 1 ，d -1 Get an order d k And a set of attitude angles with small values is used as a result of the positive kinematic solution at the current moment.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010697584.8A CN111914411B (en) | 2020-07-20 | 2020-07-20 | Full-attitude four-axis turntable frame angle instruction resolving method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010697584.8A CN111914411B (en) | 2020-07-20 | 2020-07-20 | Full-attitude four-axis turntable frame angle instruction resolving method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN111914411A CN111914411A (en) | 2020-11-10 |
| CN111914411B true CN111914411B (en) | 2022-11-11 |
Family
ID=73280651
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010697584.8A Active CN111914411B (en) | 2020-07-20 | 2020-07-20 | Full-attitude four-axis turntable frame angle instruction resolving method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111914411B (en) |
Families Citing this family (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112484720B (en) * | 2020-11-17 | 2023-04-04 | 天津津航计算技术研究所 | double-Euler full-attitude calculation method based on strapdown inertial navigation |
| CN114355787A (en) * | 2021-09-17 | 2022-04-15 | 北京星途探索科技有限公司 | Shaft-lacking turntable semi-physical simulation verification technology based on certain type supersonic cruise target |
| CN113934155A (en) * | 2021-09-17 | 2022-01-14 | 北京星途探索科技有限公司 | Semi-physical simulation method for verifying 3-2-1 rotation flight motion model by using vertical three-axis turntable |
| CN114184210A (en) * | 2021-12-03 | 2022-03-15 | 江西洪都航空工业集团有限责任公司 | Simulation test method, device and system based on horizontal turntable |
| CN114518709B (en) * | 2022-01-26 | 2024-05-17 | 哈尔滨工业大学 | Method, equipment and medium for resolving frame angle instruction of full-attitude four-axis turntable |
| CN114636357B (en) * | 2022-03-31 | 2023-11-10 | 北京中科宇航技术有限公司 | Aiming test design method aiming at shaking state of vertical turntable |
| CN116090097B (en) * | 2022-12-30 | 2024-07-09 | 北京机电工程研究所 | Near-water surface fluid-solid coupling finite element efficient calculation method based on equivalent water collision design |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2016198873A (en) * | 2015-04-14 | 2016-12-01 | トヨタ自動車株式会社 | Optimum control device, optimum control method, and optimum control program |
| CN107247157A (en) * | 2017-05-10 | 2017-10-13 | 哈尔滨工程大学 | Change the acquisition methods of Eulerian angles in a kind of quaternary number full-shape domain towards big attitude maneuver |
| CN109634293A (en) * | 2018-12-05 | 2019-04-16 | 浙江大学 | A kind of fixed-wing unmanned plane roller flowing control method |
Family Cites Families (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7680547B2 (en) * | 2006-06-08 | 2010-03-16 | Liu Hugh H T | Method, system and computer program for generic synchronized motion control for multiple dynamic systems |
| WO2016195852A1 (en) * | 2015-06-02 | 2016-12-08 | The Charles Stark Draper Laboratory, Inc. | Rapid slew and settle systems for small satellites |
| CN107685330B (en) * | 2017-10-18 | 2018-12-18 | 佛山华数机器人有限公司 | A kind of Inverse Kinematics Solution method for solving of six degree of freedom wrist bias series robot |
| CN108942942B (en) * | 2018-08-16 | 2020-01-07 | 居鹤华 | Multi-axis robot inverse kinematics modeling and resolving method based on axis invariants |
-
2020
- 2020-07-20 CN CN202010697584.8A patent/CN111914411B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2016198873A (en) * | 2015-04-14 | 2016-12-01 | トヨタ自動車株式会社 | Optimum control device, optimum control method, and optimum control program |
| CN107247157A (en) * | 2017-05-10 | 2017-10-13 | 哈尔滨工程大学 | Change the acquisition methods of Eulerian angles in a kind of quaternary number full-shape domain towards big attitude maneuver |
| CN109634293A (en) * | 2018-12-05 | 2019-04-16 | 浙江大学 | A kind of fixed-wing unmanned plane roller flowing control method |
Non-Patent Citations (5)
| Title |
|---|
| Multi-DOF Motion control system design and realization based on etherCAT;Baoxiang Xing;《IEEE》;20161208;1-7 * |
| Using dual euler angles for the analysis of arm movement during the badminton smash;Koon kiat teu;《Sports engineering》;20050131;171-178 * |
| 基于位置域迭代学习的激光导引头测试系统时变周期干扰抑制;邢宝祥;《红外与激光工程》;20190925(第2019年09期);240-247 * |
| 基于约束最优化理论的四轴转台框架角解算方法研究;徐勤贝;《中国优秀硕士学位论文全文数据库》;20170215;C032-210 * |
| 飞行器半实物仿真装备研究进展与展望;杨宝庆;《宇航学报》;20200630(第2020年06期);657-665 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN111914411A (en) | 2020-11-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN111914411B (en) | Full-attitude four-axis turntable frame angle instruction resolving method | |
| CN107688295B (en) | Four-rotor aircraft finite time self-adaptive control method based on rapid terminal sliding mode | |
| CN105159305B (en) | A kind of quadrotor flight control method based on sliding moding structure | |
| CN109062043B (en) | Spacecraft active disturbance rejection control method considering network transmission and actuator saturation | |
| CN108363298A (en) | A kind of quadrotor drone Fast Convergent control method based on quaternion representation | |
| CN106777489A (en) | UAV system opto-electric stabilization turntable tracks state modeling and simulating method | |
| CN108594837A (en) | Model-free quadrotor drone contrail tracker and method based on PD-SMC and RISE | |
| CN111695193B (en) | Modeling method and system of globally relevant three-dimensional aerodynamic mathematical model | |
| CN108845588A (en) | A kind of quadrotor Trajectory Tracking Control method based on Nonlinear Guidance | |
| CN111624878B (en) | Integral sliding mode acquisition method and system for autonomous water surface robot trajectory tracking | |
| CN112327926B (en) | Self-adaptive sliding mode control method for unmanned aerial vehicle formation | |
| CN108427428B (en) | Self-adaptive sliding mode variable structure spacecraft attitude control method based on improved iterative algorithm | |
| CN111238489B (en) | Low-earth-orbit satellite atmospheric resistance perturbation modeling and calculating method | |
| CN111880552B (en) | Multi-rotor unmanned aerial vehicle trajectory tracking composite control method | |
| CN115366109A (en) | Composite layered anti-interference method for rotor flight mechanical arm | |
| CN113253617A (en) | Online self-adaptive control method for quad-rotor unmanned aerial vehicle | |
| CN112683261A (en) | Unmanned aerial vehicle robustness navigation method based on speed prediction | |
| CN113867143A (en) | Extraterrestrial celestial body safety soft landing analysis obstacle avoidance guidance method | |
| CN117289709B (en) | High-ultrasonic-speed appearance-changing aircraft attitude control method based on deep reinforcement learning | |
| CN116360258A (en) | Hypersonic deformed aircraft anti-interference control method based on fixed time convergence | |
| CN114265420B (en) | Guidance control integrated design method suitable for high dynamic and slow response control | |
| CN112506209B (en) | Reentry vehicle prediction control method based on self-adaptive prediction period | |
| Chen et al. | Take-off and landing control for a coaxial ducted fan unmanned helicopter | |
| CN113885549A (en) | Four-rotor attitude trajectory control method based on dimension cutting PPO algorithm | |
| CN113985732A (en) | Adaptive neural network control method and device for aircraft system |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |