A method for efficiently and stably simulating non-stretchable ropes in a crane or hoist

Abstract

The present invention discloses a high-efficiency and stable method for simulating inextensible ropes in crane and hoist. The piecewise linear inextensible rope is usually represented by the Cartesian coordinates of the vertices and the quaternions on the segments, which uses too many degrees of freedom and needs many constraints, bringing unnecessary numerical difficulties and computational burden to the simulation. The present invention proposes a compact representation which uses the minimum number of degrees of freedom and naturally satisfies all the additional constraints. Specifically, the rope is regarded as a chain of rigid segments, and its shape is encoded as the Cartesian coordinates of its root vertex and the axis-angle representation of the local coordinate system on each segment. Under the representation of the present invention, the implicit time-stepping matrix has a special non-zero pattern. Using the non-zero pattern, the present invention designs a preconditioner which can solve the related linear equations with approximate linear complexity, and its speed is improved by one to two orders of magnitude compared with the widely used block-diagonal linear PCG solver.

Description

Technical Field

[0001] This invention belongs to the field of mechanical engineering, specifically relating to an efficient and stable method for simulating non-stretchable ropes in cranes and hoists. Background Technology

[0002] Hook stabilization control is a common operation in crane and hoisting work. Manual control is often very time-consuming and requires extensive operator training. Now, based on advancements in control algorithms and reinforcement learning, automatic hook stabilization technology is emerging. Among these technologies, forward simulation is the most expensive. Traditional piecewise linear non-scalable ropes are typically represented using Cartesian coordinates of vertices and quaternions on segments. This method uses too many degrees of freedom; to describe an accurate inextensible rope, three additional types of constraints need to be imposed: inextensibility constraints, unit quaternion constraints, and constraints aligning the local coordinate system axes with the rope. Each segment of the discretized rope requires these three types of constraints, resulting in a huge number of constraints. Common methods to satisfy these constraints include the penalty function method and the Lagrange multiplier method, both of which introduce unnecessary numerical difficulties and computational burdens to the simulation. Specifically, the penalty function method has three major drawbacks: it not only introduces additional numerical dissipation into the dynamic system, but also requires manual adjustment of the penalty function's weighting coefficients, and it does not accurately satisfy the aforementioned three types of constraints. The Lagrange multiplier method has two major drawbacks: it not only results in a huge system matrix, significantly increasing the solution cost, but also the fact that the three types of constraints mentioned above are all nonlinear constraints, requiring multiple iterations to accurately satisfy them. Therefore, this invention provides a highly efficient and stable method for simulating inextensible ropes. Summary of the Invention

[0003] To address the aforementioned problems, this invention proposes an efficient and stable method for simulating ropes in cranes and hoists. In this method, the rod is treated as a rigid chain of segments, and its shape is encoded by the Cartesian coordinates of its root vertex and the axis-angle representation of the material frame on each segment. Under this representation, the implicitly time-stepped Hessian matrix has a special non-zero arrangement. Utilizing this property, the relevant linear equations can be solved with approximately linear complexity. Furthermore, this invention also designs a precondition, which improves the computational speed by one or two orders of magnitude.

[0004] To achieve the above objectives, the technical solution adopted by the present invention is as follows:

[0005] An efficient and stable method for simulating non-stretchable ropes in cranes and hoists includes:

[0006] S1 represents the inextensible rope as the Cartesian coordinates of the starting point and the axis angle of each segment, thus obtaining the transformation relationship between the rope's degrees of freedom and its compact representation.

[0007] S2, model the dynamic equation of the rope based on the degree of freedom of the rope, use implicit Euler to discretize the dynamic equation of the rope in time, and obtain time discrete points with an interval of time step Δt. Each discrete point corresponds to one frame.

[0008] S3 transforms the time-discrete dynamic equations into an optimization problem, and uses a sequential quadratic programming method to solve the optimization problem, obtaining a compact expression of the rope for each frame, thus simulating an instretchable rope.

[0009] Furthermore, step S1 specifically includes:

[0010] S11, the Cartesian coordinates of the starting point of the inextensible rope are represented as p0 = {p 0，x p 0，y p 0，z Let ω represent the axis angle of the i-th segment of the inextensible rope. i ={ω i，0 ω i，1 ω i，2};

[0011] Where i = 1, 2, ..., n-1, n represents the number of vertices in the rope, p 0，x p 0，y p 0，z These are the x, y, and z coordinates of the starting point p0, respectively, and ω. i，0 ω i，1 ω i，2 These are the three components of the i-th segment's axial angle;

[0012] S12, quaternion expression for each segment of the rope calculated based on the axis angle:

[0013]

[0014] Where, q i Let ||.|| represent the quaternion representation of the i-th segment, ||.|| denotes modulo, and ω. i This represents the axial angle of the i-th segment;

[0015] S13, Calculate the local coordinate system based on the axis angle:

[0016] m i =exp ω i

[0017] Where, m i Let i represent the local coordinate system of the i-th segment;

[0018] S14, based on the Cartesian coordinates of the starting point, recursively deduce the Cartesian coordinates of the remaining vertices. The recursive calculation formula is:

[0019] pi+1 =p i +l i d(m i )

[0020] Where, p i Let l represent the Cartesian coordinates of the (i-1)th vertex. i Let d represent the length of the i-th segment, and d(.) represent the extraction of the third axis in the local coordinate system.

[0021] S15, the degrees of freedom of the rope are expressed as:

[0022]

[0023]

[0024] p0={p 0，x p 0，y p 0，z}

[0025] ω=(ω1, ω2, ..., ω i ,...,ω n-1 )

[0026] ω i ={ω i，0 ω i，1 ω i，2}

[0027] Where x represents the degree of freedom of the rope, p represents the position of each vertex in the rope, q represents the quaternion representation of each segment in the rope; T(.) represents the coordinate transformation function, used to represent the transformation relationship between the degree of freedom of the rope and the compact representation of the rope; s represents the compact representation of the rope, ω represents the set of axis angles of each segment in the rope, and the superscript T represents transpose.

[0028] Furthermore, step S3 specifically includes:

[0029] S31, the time-discrete dynamic equations are transformed into objective functions of bending torsional potential energy V(q) and translational kinetic energy increment T. p (p), rotational kinetic energy T q The sum of (q) and gravitational potential energy G(p), and subject to an impenetrable constraint C from the collision. I The optimization problem with (x)≥0 is represented as:

[0030] x t+Δt =argminE(x)=V(q)+T p (p)+T q (q)+G(p)

[0031] stC I(x)≥0

[0032] Where, x t+Δt Let E(x) represent the degrees of freedom of the rope in the next frame, and let C represent the objective function after implicit Euler transformation. I (x) represents the distance between points, planes, or segments, x represents the degree of freedom of the rope, p represents the position of each vertex in the rope, and q represents the quaternion representation of each segment in the rope;

[0033] S32, transform the optimization problem in step S31 into an optimization problem where the degrees of freedom x are replaced by a compact representation s of the rope:

[0034] s t+Δt =argminE(T(s))

[0035] stC I (T(s))≥0

[0036] Where s represents the compact representation of the rope, and T(.) represents the coordinate transformation function, used to represent the transformation relationship between the rope's degrees of freedom and its compact representation. t+Δt This indicates a compact representation of the rope in the next frame;

[0037] S33, The line search method is selected to ensure convergence. The metric function φ(.) for the line search method is:

[0038]

[0039] Where μ represents the reciprocal of the maximum value of the Lagrange multiplier, and ||.|| represents the modulo operation;

[0040] S34, Solve the optimization problem using sequential quadratic programming:

[0041]

[0042]

[0043]

[0044]

[0045] Where, Δs k Let A represent the increment of the compact expression s in the k-th iteration, where the superscript T indicates transpose, Δs represents the increment of the compact expression s, and A represents the increment of the compact expression s. k This represents the system matrix in the k-th iteration. This indicates that the derivative of the compact expression s, b k Let J(.) denote the right-hand term of the k-th iteration, and let J(.) denote the derivative of the coordinate transformation function T(.). kLet H represent the compact expression of the k-th iteration, and let H represent the Hessian matrix of the objective function E(x) with x as the variable.

[0046] Furthermore, the bending torsional potential energy V(q) and the translational kinetic energy increment T... p (p), rotational kinetic energy T q The formulas for calculating gravitational potential energy G(p) and gravitational potential energy G(q) are as follows:

[0047]

[0048]

[0049]

[0050]

[0051] Among them, M p For the lumped mass matrix, p t Let p be the Cartesian coordinate of each vertex in the previous frame, and Δt be the time step. J is the velocity of each vertex in the previous frame. kk For the moment of inertia, l i Let B be the length of the i-th segment. k Let q be a constant antisymmetric matrix. i Represent the quaternion representation of the i-th segment. Let h represent the time derivative of the quaternion expression for the i-th segment. i Let p be the mass of the i-th vertex, g be the acceleration due to gravity, and p be the mass of the i-th vertex. i，y Let K be the y-coordinate of the i-th vertex. kk This is the diagonal matrix of element stiffness.

[0052] Furthermore, in step S34, when solving the optimization problem using the sequential quadratic programming method, the active set method is used to handle collisions and generate the corresponding KKT matrix, and the Schur complement theory method is used for solving, specifically:

[0053] First, collision pairs are obtained using the discrete collision detection method, and point-to-surface collision constraints and segment-to-segment collision constraints C are generated. I The expression is as follows:

[0054]

[0055] Where, p a Let p be the coordinates of the point in the point-to-surface collision. i p j and p k Let p be the coordinates of three vertices of a triangular facet in a point-to-face collision. A segment-to-segment collision pair consists of segments from two ropes, where one segment has two vertices p.i and p i+1 The two vertices of the other segment are p. j and p j+1 w i w j w k Represented as p i p j p k The centroid coordinates, δ represents the collision detection threshold, and N represents the collision normal;

[0056] Then, the active set method is used to extract C. I Active collision constraints C in A , The quadratic programming problem in step S34 is transformed into solving the following KKT matrix:

[0057]

[0058] Where λ represents the Lagrange multiplier;

[0059] Using the Schur complement theory, the KKT matrix is ​​transformed into the following system of equations:

[0060]

[0061] Among them, (A) k ) -1 A represents k The reverse.

[0062] Furthermore, the preconditional conjugate gradient method is used to apply (A) k ) -1 To solve the problem, a tridiagonal matrix is ​​used as a precondition for A. k The equation is solved as a left-hand term, with the following preconditions:

[0063]

[0064] Where P is the precondition, and diagblk(.) is the function to extract the diagonal portion of the block; J p，s J is the derivative of the position p of each vertex in the rope with respect to the compact expression s of the rope. q，ω Let H be a quaternion expression for the derivative of q with respect to the axial angle ω of each segment of the rope. p H is the derivative of the objective function E(x) with respect to the positions p of each vertex in the rope. q Let q be the derivative of the objective function E(x) with respect to the quaternion expression of each segment of the rope, with superscript q. This indicates transpose.

[0065] Furthermore, the preconditions mentioned above and Both can be constructed in linear time complexity, and the construction method is as follows:

[0066] (1) Given J q，ω and H q In this case, The calculation method for the (i, j)th block is as follows:

[0067]

[0068] Among them, blk i，i (.), blk i，j (.), blk j，j (.) represent extracting the (i, i), (i, j), and (j, j)th blocks from the matrix, respectively.

[0069] (2) Given J p，s and H p In this case, The calculation method for the (i, i)th block is as follows:

[0070]

[0071] Among them, blk k，k (.) indicates the extraction of the (k, k)th block in the matrix, where n represents the number of vertices of the rope.

[0072] Compared with the prior art, the advantages of the present invention are:

[0073] (1) Because the present invention adopts a compact form of expression, compared with the traditional discretization method, it avoids the introduction of non-stretchable constraints, unit quaternion constraints and the constraints of the alignment of the axes and ropes of the local coordinate system, thereby greatly reducing the condition number of the dynamic equation.

[0074] (2) Because the present invention adopts a compact expression form, compared with the penalty function method under the traditional discrete method, it avoids the introduction of a large weight penalty term, reduces the numerical dissipation of the simulation, and can accurately characterize the inextensible rope.

[0075] (3) Because the present invention adopts a compact expression form, compared with the traditional discrete method of Lagrange multipliers, it avoids introducing too many constraint terms, thereby greatly simplifying the final KKT matrix and improving the efficiency of solving.

[0076] (4) Because this invention explores the special structure of the objective function Hessian matrix under this compact expression and develops a precondition with a tridiagonal structure that can be constructed with linear complexity, compared with the existing diagonal preconditions, it better preserves the terms of the original matrix, thereby improving the computational complexity by one to two orders of magnitude. Attached Figure Description

[0077] Figure 1 This is a schematic diagram illustrating the discrete representation of a rope in one embodiment of the present invention.

[0078] Figure 2 This is a flowchart of the sequential quadratic programming algorithm in one embodiment of the present invention.

[0079] Figure 3 This is a schematic diagram of the calculation results obtained in one embodiment of the present invention.

[0080] Figure 4 This is a schematic diagram of the calculated results in another embodiment of the present invention. Detailed Implementation

[0081] To make the above-mentioned objects, features, and advantages of the present invention more apparent and understandable, the specific embodiments of the present invention will be described in detail below with reference to the accompanying drawings. Many specific details are set forth in the following description to provide a thorough understanding of the present invention. However, the present invention can be practiced in many other ways different from those described herein, and those skilled in the art can make similar modifications without departing from the spirit of the present invention. Therefore, the present invention is not limited to the specific embodiments disclosed below. Technical features in the various embodiments of the present invention can be combined accordingly without mutual conflict.

[0082] In the description of this invention, it should be understood that the terms "first" and "second" are used only for descriptive purposes and should not be construed as indicating or implying relative importance or implicitly specifying the number of indicated technical features. Therefore, a feature defined with "first" and "second" may explicitly or implicitly include at least one of those features.

[0083] The present invention proposes a highly efficient and stable method for simulating inextensible ropes, which mainly includes the following steps:

[0084] Step 1: Discretize the instretchable rope as the Cartesian coordinates of the starting point and the axis angle of each segment.

[0085] like Figure 1 As shown, for a given array of n vertices (p0, p1, ..., p...), ... n-2 p n-1 The discrete rope has a total of n-1 segments. The compact representation proposed in this invention uses only the Cartesian coordinates of the first point p0 and the axis-angle representation of each segment. The direction of each segment is obtained through the axis-angle representation, thereby recursively obtaining the direction of the second vertex p1 to the nth vertex p. n-1 Cartesian coordinates.

[0086] S11. Represent the starting point as p. c =p0={p0，x p 0，y p 0，z Let the axial angle of the i-th segment be ω. i ={ω i，0 ω i，1 ω i，2}, i = 1, 2, ..., n-1, p 0，x p 0，y p 0，z These are the x, y, and z coordinates of point p0, respectively, and ω. i，0 ω i，1 ω i，2 These are the three components of the axial angle of the i-th segment.

[0087] S12. Calculate the quaternion expression for each segment based on the axis angle:

[0088]

[0089] Where, q i Let ||.|| represent the quaternion representation of the i-th segment, ||.|| denotes modulo, and ω. i This represents the axial angle of the i-th segment;

[0090] S13. Calculate the local coordinate system based on the axis angle;

[0091] m i =exp ω i

[0092] Where, m i Let i represent the local coordinate system of the i-th segment;

[0093] S14. Based on the Cartesian coordinates of the starting point, recursively deduce the Cartesian coordinates of the remaining vertices. The recursive calculation formula is as follows:

[0094] p i+1 =p i +l i d(m i )

[0095] Where, p i Let l represent the Cartesian coordinates of the (i-1)th vertex. i Let d represent the length of the i-th segment, and d(.) represent the extraction of the third axis in the local coordinate system.

[0096] S15. Express the degrees of freedom of the rope as:

[0097]

[0098]

[0099] p0={p 0，x p0，y p 0，z}

[0100] ω=(ω1, ω2, ..., ω i ,...,ω n-1 )

[0101] ω i ={ω i，0 ω i，1 ω i，2}

[0102] Where x represents the degrees of freedom of the rope, p represents the position of each vertex in the rope, q represents the quaternion representation of each segment in the rope, T(.) represents the coordinate transformation function, s represents the compact expression of the rope, ω represents the set of axes and angles of each segment in the rope, and the superscript T represents the transpose.

[0103] Step 2: Model the dynamic equation of the rope based on the degree of freedom of the rope. Use implicit Euler to discretize the dynamic equation of the rope in time, and obtain time discrete points with an interval of time step Δt. Each discrete point corresponds to one frame. Transform the time-discrete dynamic equation into an optimization problem. Use the sequential quadratic programming method to solve the optimization problem and obtain a compact expression of the rope corresponding to each frame, thus realizing the simulation of an inextensible rope.

[0104] The specific implementation method in this step is as follows:

[0105] S21. Modeling based on Kirchhoff's theory The dynamic equations of a rope with degrees of freedom are transformed into objective functions: bending and torsional potential energy V(q) and translational kinetic energy increment T. p (p), rotational kinetic energy T q The sum of (q) and gravitational potential energy G(p) and subject to an impenetrable constraint C from the collision. I Optimization problem where (x)≥0:

[0106] x t+Δt =argminE(x)=V(q)+T p (p)+T q (q)+G(p)

[0107] stC I (x)≥0

[0108] Where, x t+Δt Let E(x) represent the degrees of freedom of the rope in the next frame, and let C represent the objective function after implicit Euler transformation. I (x) represents the distance between points, surfaces, or segments, x represents the degree of freedom of the rope, p represents the position of each vertex in the rope, and q represents the quaternion representation of each segment in the rope.

[0109] The formula for calculating the increase in translational kinetic energy is as follows:

[0110]

[0111] Among them, M p For the lumped mass matrix, p t Let p be the Cartesian coordinate of each vertex in the previous frame, and Δt be the time step. This represents the velocity of each vertex in the previous frame.

[0112] The formula for calculating rotational kinetic energy is:

[0113]

[0114] Among them, f kk For rotational inertia, The time derivative of a quaternion, l i Let B be the length of the i-th segment. k Let q be a constant antisymmetric matrix. i Describe the quaternion of the i-th segment. Let represent the time derivative of the quaternion in the i-th segment.

[0115] The formula for calculating bending torsional potential energy is:

[0116]

[0117] in, K represents the natural curvature of the rope. kk Let K be the element stiffness diagonal matrix. The terms of this matrix are: 11 =K 22 =Eπr 2 / 4,K 33 =Gπr 2 / 2. Where E and G are Young's modulus and shear modulus, respectively, and r is the cross-sectional radius of the rope; B k It is a constant antisymmetric matrix, with the following specific form:

[0118]

[0119] The formula for calculating gravitational potential energy is:

[0120]

[0121] Among them, h i Let p be the mass of the i-th vertex, g be the acceleration due to gravity (assuming gravity is applied in the direction of g), and p be the acceleration due to gravity. i，y Let y be the y-coordinate of the i-th vertex.

[0122] S22, transform the optimization problem in step S21 into an optimization problem where the degrees of freedom x are replaced by a compact representation s of the rope:

[0123] s t+Δt =argminE(T(s))

[0124] stC I (T(s))≥0

[0125] Where s represents the compact representation of the rope, and T(.) represents the coordinate transformation function, used to represent the transformation relationship between the rope's degrees of freedom and its compact representation. t+Δt This indicates a compact representation of the rope in the next frame;

[0126] S23, The line search method is selected to ensure convergence. The metric function φ(.) for the line search method is:

[0127]

[0128] Where μ represents the reciprocal of the maximum value of the Lagrange multiplier, and ||.|| represents the modulo operation; in this embodiment, the search step size α is calculated using the line search method;

[0129] S24, Solve the optimization problem using sequential quadratic programming:

[0130]

[0131]

[0132] Where, Δs k Let A represent the increment of the compact expression s in the k-th iteration, where the superscript T indicates transpose, Δs represents the increment of the compact expression s, and A represents the increment of the compact expression s. k This represents the system matrix in the k-th iteration. Let denote the derivative of the compact expression s, where bk represents the right-hand term of the k-th iteration; the system matrix A k and right-hand item b k The calculation formula is:

[0133]

[0134]

[0135] Where J(.) represents the derivative of the coordinate transformation function T(.), s k Let H represent the compact expression of the k-th iteration, and let H represent the Hessian matrix of the objective function E(x) with x as the variable.

[0136] In a specific embodiment of the present invention, when solving the quadratic programming problem described in step S24 above, the active set method is used to handle collisions and generate the corresponding KKT matrix, and the Schur Complement method is used for solving.

[0137] The specific implementation steps are as follows:

[0138] S241, collision pairs are obtained according to the discrete collision detection method, and point-to-surface collision constraints and segment-to-segment collision constraints C are generated. I as follows:

[0139]

[0140] Where, p a Let p be the coordinates of the point in the point-to-surface collision. i p j and p k Let be the coordinates of the three vertices of the triangular facet in a point-to-face collision.

[0141] A segment collision pair consists of segments in two ropes, where one segment has two vertices p. i and p i+1 The two vertices of the other segment are p. j and p j+1 w i w j w k Represented as p i p j p k The centroid coordinates, δ represents the collision detection threshold, and N represents the collision normal;

[0142] S242, using the active set method to extract C I Active collision constraints C in A ,

[0143] S243, the quadratic programming problem in S24 is transformed into solving the following KKT matrix:

[0144]

[0145] in, For active collision constraints, λ denotes the Lagrange multiplier;

[0146] S244, this matrix is ​​solved using the Schur Complement method, which is equivalent to solving the following two systems of equations:

[0147]

[0148] Where λ represents the Lagrange multiplier, (A k ) -1 A represents k The inverse of A is not directly calculated in the implementation, but rather A is... k Therefore, when solving this system of equations, the most expensive step is to solve the equation with A as the left-hand term. k The equations for the left-hand terms require solving a total of +2 equations with constraints.

[0149] To accelerate the solution, this invention uses the preconditional conjugate gradient method (PCG). Specifically, for A... k Using a special tridiagonal matrix as a precondition, the speed is accelerated to A k The solution to the equation for the left-hand term is as follows:

[0150] First, according to step S24, A = J T HJ, matrix H has the following structure:

[0151]

[0152]

[0153]

[0154] Among them, H p H is the derivative of the objective function E(x) with respect to the positions p of each vertex in the rope. q Let K be the derivative of the objective function E(x) with respect to the quaternion expression q of each segment of the rope. q M represents the Hessian matrix of the bending torsional potential energy with respect to q. q This represents the generalized rotational inertia matrix.

[0155] Matrix J has the following structure:

[0156]

[0157] Among them, J q，s Let q be the derivative of the quaternion expression q of each segment of the rope with respect to the compact expression s of the rope. This matrix is ​​a block diagonal matrix, specifically representing the sound... =0, J is a block diagonal matrix; p，s Let p be the derivative of the position p of each vertex in the rope with respect to the compact representation s of the rope. This matrix is ​​a lower triangular matrix.

[0158] Subsequently, based on the above special structure, the system matrix A can be simplified to:

[0159]

[0160] The precondition is:

[0161]

[0162] Where P is a precondition. To extract The diagonal part of the block, The quaternion representing each segment of the rope expresses the derivative of q with respect to the axial angle ω of each segment of the rope;

[0163] The above preconditions and Both can be constructed with linear complexity, therefore P can be constructed with linear complexity.

[0164] In one specific embodiment of the invention, a linear construction and The method is as follows:

[0165] Given J q，ω and H q In this case, The calculation method for the (i, j)th block is as follows:

[0166]

[0167] Among them, blk i，i (.), blk i，j (.), blk j，j (.) represent extracting the (i, i), (i, j), and (j, j)-th blocks from the matrix, respectively; because The result is a tridiagonal matrix because the overall construction complexity is linear.

[0168] Given J p，s and H p In this case, The calculation method for the (i, i)th block is as follows:

[0169]

[0170] Among them, blk k，k (.) indicates extracting the (k, k)th block from the matrix, where n represents the number of vertices in the rope. The recursive calculation formula achieves linear complexity.

[0171] Figure 2 This is a flowchart of the simulation algorithm of the present invention, which uses a given time step Δt and the current state x. k Calculate the state x after a step size of Δt. k+1 The specific process is as follows:

[0172] Step 101, calculate the degrees of freedom x of the rope in the k-th iteration. k As the state of the previous frame, x k The time step Δt is used as input.

[0173] Step 102, based on the state x of the previous frame k Calculate the compact representation of the rope in the previous frame. k :

[0174] s k =T -1 (x k )

[0175] Step 103, determine whether the current state has reached instantaneous force equilibrium:

[0176]

[0177] If a force balance is achieved, proceed to step 108; otherwise, proceed to step 104.

[0178] Step 104: Calculate the system matrix A for the quadratic programming problem. k Right-hand item b k and constraint matrix

[0179] Step 105: Use the active set method to solve for the increment Δs and Lagrange multiplier λ of the compact representation.

[0180] Step 106: Calculate the step size α using the line search method.

[0181] Step 107, calculate the compact representation s of the rope in the next frame. k+1 =s k +αΔs.

[0182] Step 108, calculate the degree of freedom x of the rope in the next frame. k+1 =T(s) k+1 (This is used as the state for the next frame.)

[0183] Figure 3 This is a schematic diagram illustrating the result of calculating a rope knot using the method of this invention. Figure 3 From left to right, the frames represent frames 0 (initial state), 40, 120, 195, and 273 of the simulation, with a time step of 0.003 seconds. The rope is discretized into 101 points and 100 segments. The Young's modulus is 20 MPa, the shear modulus is 20 MPa, the cross-sectional radius of the rope is 0.03 m, and the density is 0.1 kg / m³. 3In this embodiment, a tension of 200N is applied to both ends of the rope, forming a knot. This invention achieved a frame rate of 25fps on a desktop computer with an AMD Ryzen 7-3800X 3.9GHz CPU and 32GB of RAM.

[0184] Figure 4 The method of this invention is used to calculate the process of twisting a straight, inextensible rope into a double helix structure. Figure 4 From left to right, the frames are 200, 260, 700, and 1999, with a simulation time step of 0.01 seconds. The rope is discretized into 51 points and 50 segments. The Young's law of thermal expansion is 1000 MPa, the shear modulus is 100000 MPa, the cross-sectional radius is 0.003 m, and the density is 200 kg / m³. 3 By applying two opposing torques to both ends of the rope and restricting movement to the axial direction only, the rope slowly forms a spiral structure. This invention achieved a frame rate of 50fps on a desktop computer with an AMD Ryzen 7-3800X 3.9GHz CPU and 32GB of RAM.

[0185] The above examples are merely specific embodiments of the present invention. Obviously, the present invention is not limited to the above embodiments and many variations are possible. All variations that can be directly derived or conceived by those skilled in the art from the disclosure of the present invention should be considered within the scope of protection of the present invention.

Claims

1. A highly efficient and stable method for simulating non-stretchable ropes in cranes and hoists, characterized in that, include: S1 represents the inextensible rope as the Cartesian coordinates of the starting point and the axis angle of each segment, thus obtaining the transformation relationship between the rope's degrees of freedom and its compact representation. S2, Model the dynamic equations of the rope based on the degree of freedom representation of the rope, and use implicit Euler discretization to discretize the dynamic equations of the rope in time, obtaining the interval as the time step. The time discrete points, each discrete point corresponds to one frame; S3 transforms the time-discrete dynamic equations into an optimization problem, and uses a sequential quadratic programming method to solve the optimization problem, obtaining a compact expression of the rope for each frame, thus simulating an instretchable rope.

2. The method for efficiently and stably simulating non-stretchable ropes in cranes and hoists according to claim 1, characterized in that, The specific steps of S1 are as follows: S11, representing the Cartesian coordinates of the starting point of the instretchable rope as follows: The first of the non-stretchable rope The axial angle of the segment is expressed as ; in, , Indicates the number of vertices in the rope. These are the starting points. x, y, z coordinates They are the first The three components of the segment axis angle; S12, quaternion expression for each segment of the rope calculated based on the axis angle: in, Indicates the first Quaternion representation of segments Indicates modulo, Indicates the first The segment's axial angle is represented; S13, Calculate the local coordinate system based on the axis angle: in, Indicates the first The local coordinate system of the segment; S14, based on the Cartesian coordinates of the starting point, recursively deduce the Cartesian coordinates of the remaining vertices. The recursive calculation formula is: in, Indicates the first Cartesian coordinates of the vertices, Indicates the first The length of the segment This indicates that the third axis in the local coordinate system is extracted; S15, the degrees of freedom of the rope are expressed as: Where x represents the degree of freedom of the rope, This indicates the position of each vertex in the rope. The quaternion representation of each segment of the rope; This represents a coordinate transformation function used to express the transformation relationship between the degrees of freedom of a rope and its compact representation; A compact expression indicating the rope. Represents the set of axial angles of each segment of a rope, with superscript. This indicates transpose.

3. The method for efficiently and stably simulating non-stretchable ropes in cranes and hoists according to claim 1, characterized in that, The specific steps of step S3 are as follows: S31, transforming the time-discrete dynamic equations into an objective function of bending torsional potential energy. Translational kinetic energy increment Rotational kinetic energy and gravitational potential energy The sum, and subject to impenetrable constraints from collisions. The optimization problem is represented as: in, This indicates the degree of freedom of the rope in the next frame. Let represent the optimization objective function after implicit Euler transformation. This represents the distance between points, planes, or segments; x represents the degrees of freedom of the rope. This indicates the position of each vertex in the rope. The quaternion representation of each segment of the rope; S32 transforms the optimization problem in step S31 into a compact expression using ropes. The optimization problem of replacing the degree of freedom x: in, A compact expression indicating the rope. This represents a coordinate transformation function used to express the transformation relationship between the degrees of freedom of a rope and its compact representation. This indicates a compact representation of the rope in the next frame; S33, The line search method is selected to ensure convergence, and the metric function of the line search method is... for: Where µ represents the reciprocal of the maximum value of the Lagrange multiplier. Indicates modulo; S34, Solve the optimization problem using sequential quadratic programming: in, Indicates the first Compact expression of the next iteration The increment, superscript Indicates transpose. Indicates a compact expression The increment, Indicates the first The system matrix of the next iteration Indicates a compact expression Differentiate, Indicates the first The right-hand term of the next iteration Represents coordinate transformation function The derivative of This represents a compact representation of the k-th iteration. Describe the objective function by The Hessian matrix with variables.

4. The method for efficiently and stably simulating non-stretchable ropes in cranes and hoists according to claim 3, characterized in that, The bending torsional potential energy Translational kinetic energy increment Rotational kinetic energy and gravitational potential energy The calculation formula is: in, Indicates the natural curvature of the rope. For lumped mass matrix, Let be the Cartesian coordinates of each vertex in the previous frame. Let Cartesian coordinates be the coordinates of each vertex. For time step, The velocity of each vertex in the previous frame. For rotational inertia, For the first The length of the segment It is a constant antisymmetric matrix. Indicates the first Quaternion representation of segments Indicates the first The time derivative of the segment expressed by quaternions. For the first The mass of each vertex It is the acceleration due to gravity. For the first vertices coordinate, This is the diagonal matrix of element stiffness.

5. The method for efficiently and stably simulating non-stretchable ropes in cranes and hoists according to claim 3, characterized in that, In step S34, when solving the optimization problem using the sequential quadratic programming method, the active set method is used to handle collisions and generate the corresponding KKT matrix, and the Schur complement theory method is used for solving, specifically: First, collision pairs are obtained using the discrete collision detection method, and point-to-surface collision constraints and segment-to-segment collision constraints are generated. The expression is as follows: in, Let these be the coordinates of the point in the point-to-surface collision. , and Let G be the coordinates of three vertices of a triangular facet in a point-to-face collision. In a segment-to-segment collision pair, there are segments from two ropes, where one segment has two vertices... and The two vertices of the other segment are and , , They are respectively represented as , , The coordinates of the center of gravity, This represents the threshold for collision detection. Indicates the collision normal; Then, the active set method is used to extract... Active collision constraints in , The quadratic programming problem in step S34 is transformed into solving the following KKT matrix: in, Represents the Lagrange multipliers; Using the Schur complement theory, the KKT matrix is ​​transformed into the following system of equations: in, express The reverse.

6. The method for efficiently and stably simulating non-stretchable ropes in cranes and hoists according to claim 5, characterized in that, Using the preconditional conjugate gradient method To solve the problem, a tridiagonal matrix is ​​used as a precondition during the solution process. The equation is solved as a left-hand term, with the following preconditions: in, As a prerequisite, A function to extract the diagonal portion of a block; The positions of the vertices in the rope A compact expression of rope The derivative of Quaternion representation of each segment of the rope The axial angle of each segment of the rope The derivative of For the objective function Position of each vertex in the rope The derivative of For the objective function Quaternion representation of each segment of the rope The derivative, superscript This indicates transpose.

7. The method for efficiently and stably simulating non-stretchable ropes in cranes and hoists according to claim 6, characterized in that, The aforementioned preconditions and Both can be constructed in linear time complexity, and the construction method is as follows: (1) In the known and In this case, The The calculation method for each block is as follows: in, They represent the extraction matrix of the first... , , Each block; (2) In the known and In this case, The The calculation method for each block is as follows: in, Indicates the extraction matrix of the first... Each block, This indicates the number of vertices in the rope.

Images