A method for fast simulation of three-dimensional dynamic contact behavior of a pantograph-rigid catenary system

Abstract

This invention discloses a rapid simulation method for the three-dimensional dynamic contact behavior of a pantograph-rigid contact network. Specifically, it involves: establishing a rigid contact network model using beam elements and three-dimensional spring elements; establishing a pantograph reduced mass block model; and constructing the pantograph-rigid contact network motion equations. The method uses Newton's iteration method to calculate the initial equilibrium state of the system under gravity, and calculates the elastic internal force vector, tangent stiffness matrix, and generalized stiffness matrix of the system at this time. It then calculates the pantograph-contact network contact force vector at the current time step and the system coordinate vector for the next time step. Finally, it determines whether the calculations for all time steps have been completed. If completed, the contact force sequence is output and the calculation ends; otherwise, the time step index and system state are updated, and the process returns to the previous step to perform calculations for subsequent time steps. This invention significantly improves the simulation efficiency of the pantograph-rigid contact network system, significantly reduces simulation time costs, and ensures extremely high computational accuracy.

Description

Technical Field

[0001] This invention belongs to the field of dynamic simulation technology for electrified railways, and particularly relates to a rapid simulation method for the three-dimensional dynamic contact behavior of a pantograph rigid contact network. Background Technology

[0002] The pantograph-rigid contact wire system is an important way for electrified trains in tunnels to obtain electrical energy. Due to the high cost and harsh test conditions of actual line tests, the research work in the field of pantograph-contact wire dynamics is mainly based on the finite element model. Zhou Ning et al. (Zhou Ning, Zou Huan, Zou Dong, et al. Adaptability study of simulation model of urban rail transit pantograph-contact wire system[J]. Journal of Southwest Jiaotong University, 2017, 52(2): 408-423) compared the variable stiffness model and the finite element model of the rigid contact wire and proved that the results obtained by the finite element model are more accurate under high-speed conditions. C. Vera et al. (C. Vera, B. Suarez, J. Paulin, et al. Simulation model for the study of overhead rail current collector systems dynamics, focused on the design of a new conductor rail[J]. Veh. Syst. Dyn. 2006, 44: 595-614) established a finite element model of the rigid contact wire using the finite element software ANSYS. L. Chen et al. (L. Chen, F. Duan, Y. Song, et al. Assessment of dynamic interaction performance of high-speed pantograph and overhead conductor rail system[J]. IEEE Trans. Instrum. Meas. 2022, 71: 1-14) established a finite element model of a rigid contact network using the Absolute Nodal Coordinates (ANCF) method and solved the dynamic contact force of the pantograph and catenary using the Newmark method. In existing modeling and simulation studies of pantograph-rigid contact network systems, the penalty function method is used to simulate the pantograph-catenary contact behavior. The contact stiffness introduced by this method is an empirical parameter, and its value significantly affects the simulation results and is difficult to determine a suitable value. Furthermore, the penalty function method also suffers from the problem of virtual penetration between the pantograph and catenary, leading to inaccurate displacement results in the simulation. In addition, pantograph-rigid contact network simulations only consider the vertical contact force and do not consider the variation of the contact force direction in three-dimensional space.To address the above issues, L. Chen et al. (L. Chen, F. Duan, Y. Song, et al. Three-dimensional contact formulation for assessment of dynamic interaction of pantograph and overhead conductor railsystem[J]. Veh. Syst. Dyn. 2022: 1-24.) established an accurate contact model between the pantograph and rigid contact net based on the Lagrange multiplier constraint method, and introduced a three-dimensional vector definition format for the contact force in the global coordinate system. They also calculated the normal and tangential contact forces, fully considering the changes in the contact force direction caused by pantograph-net vibration. However, regardless of whether the penalty function method or the accurate contact model is used to simulate pantograph-net contact, iterative steps are required, and each iteration step requires a large number of matrix operations, resulting in long computation time. In the accurate contact model, iterative calculation of the tangential force direction is necessary, and the overall solution process of the accurate contact model is explicit, requiring a small time step to maintain stable convergence, resulting in low overall computational efficiency. The enormous computational time cost during optimization design research or batch simulation limits the development of related research. Summary of the Invention

[0003] To improve the simulation efficiency of pantograph-rigid contact network systems when using precise three-dimensional calculation models, this invention provides a rapid simulation method for the three-dimensional dynamic contact behavior of pantograph-rigid contact networks.

[0004] The present invention provides a rapid simulation method for the three-dimensional dynamic contact behavior of a pantograph rigid contact network, comprising the following steps:

[0005] Step 1: A finite element model of the rigid contact network is established using beam elements to simulate the busbar and contact wire, and spring elements to simulate the suspension structure and clamps. The pantograph model is then established using the reduced mass block method, and the motion equations for the pantograph-rigid contact network are constructed. Specifically:

[0006] S11: Based on the design parameters of the rigid contact network, the rigid contact network is discretized using beam elements (this invention is applicable to all types of beam elements) to determine the coordinate vectors of each element and node; based on beam element theory, the mass matrix, stiffness matrix, damping matrix, and load vector of each element are calculated.

[0007] S12: The suspension structure and clamps are considered as a whole and simulated using a three-dimensional spring element. The mass matrix M of the three-dimensional spring element is... s and stiffness matrix K s They are respectively:

[0008]

[0009]

[0010] Where, m s k is the equivalent mass of the suspension structure. x k y k z These are the equivalent stiffnesses of the three-dimensional spring on the three coordinate axes X, Y, and Z of the global coordinate system.

[0011] S13: Using the reduced mass block method, a mass block model of the pantograph is established. The motion equations of the pantograph three-mass block model are:

[0012]

[0013] Where m1, m2, and m3 are the masses of the bow head, upper frame, and lower frame, respectively; c1, c2, and c3 are the damping of the bow head, upper frame, and lower frame, respectively; k1, k2, and k3 are the stiffness of the bow head, upper frame, and lower frame, respectively; y1, y2, and y3 are the displacements of the bow head, upper frame, and lower frame, respectively; F L This refers to the static lifting force experienced by the pantograph.

[0014] S14: The motion equations of the pantograph-rigid contact network system are constructed using the standard finite element method as follows:

[0015]

[0016] Where M and C are the overall mass matrix and damping matrix of the pantograph-rigid contact network system, respectively; Let Q be the velocity and acceleration vectors of the system as a whole; let Q be the elastic internal force vector of the system as a whole; let P be the external load vector acting on the system; the formula for calculating the elastic internal force vector of the system is:

[0017] Q = Ke (5)

[0018] Where K is the system stiffness matrix and e is the system coordinate vector.

[0019] Step 2: Substitute the gravity load into the external load vector in the pantograph-rigid contact wire motion equation, and use Newton's iteration method to calculate the coordinate vector e0 of the system in the initial equilibrium state under gravity; based on the initial equilibrium state of the pantograph-contact wire system under gravity, pre-calculate the elastic internal force vector Q0 of the system in this state, and assemble the tangential stiffness matrix K of the system in this state. t The generalized stiffness matrix J is then calculated and sparsified. The expression for the system's generalized stiffness matrix J is:

[0020]

[0021] Here, Δt represents the iteration time step in the dynamic simulation process. The zero vector serves as a placeholder; its value needs to be updated in subsequent dynamic simulations. Memory is allocated in advance here to reduce the time overhead of memory reallocation.

[0022] Step 3: Calculate the contact force vector of the pantograph-rigid contact network system at the current time step, and the system coordinate vector at the next time step.

[0023] S31: Calculate the system's elastic internal force Q at the current time step using the concept of linearization. i And calculate the generalized load vector. The calculation formulas are as follows:

[0024] Q i =Q0+K t (e i -e0) (7)

[0025]

[0026] Where the subscript i is the time step index. At the initial time, i.e., when i = 0, e i-1 Calculated by the following formula:

[0027]

[0028] S32: Calculate the contact force direction vector n using the estimated relative velocity. i The calculation formula is:

[0029]

[0030] Where μ is the Coulomb friction coefficient between the bow and the fire net, and n n,i The normal contact force direction vector is...

[0031]

[0032] Among them, A p,i and A c,i These are the tangential vectors of the pantograph head and the contact wire, respectively.

[0033] n t,i The estimated tangential force direction vector.

[0034]

[0035] Among them, v i The estimated relative velocity between the bow head and the contact line at the contact point is given by the following formula:

[0036]

[0037] Among them, vt Let G be the vehicle's velocity vector. Let G be the contact constraint matrix, written as...

[0038]

[0039] Define the pantograph's forward direction as "forward" and the opposite direction as "rear," where n1 is the number of degrees of freedom of all nodes behind the beam elements in contact on the rigid contact network, and n2 is the number of degrees of freedom of all nodes in front of the beam elements in contact on the rigid contact network; S c S is the shape function matrix of a beam element currently in contact. p This is an auxiliary matrix used to calculate the vertical position of the bow head, written as:

[0040]

[0041] S33: Assuming the pantograph-catenary system is in contact, update the last row and last column of the system's generalized stiffness matrix and perform LU decomposition on it. The update format is as follows:

[0042]

[0043] Where Δt is the iteration time step, i is the iteration index, and T = [0 1 0] is the auxiliary vector. Performing LU decomposition on J yields the upper triangular matrix J. U and lower triangular matrix J L Satisfying J L ·J U =J.

[0044] S34: Calculate the contact force at the current time step and the system coordinate vector at the next time step. The calculation formula is as follows:

[0045]

[0046] Among them, U i It includes the system's coordinate vector and the magnitude of the contact force, that is:

[0047]

[0048] f i The total contact force is the magnitude of the total contact force at the current time step. The total contact force vector is calculated by the following formula:

[0049] F i =f i ·n i (19)

[0050] The actual contact state between the pantograph and the rigid contact wire is determined by the sign of the contact force. If f i If f ≥ 0, then the bow and fire net are in contact, proceed to step 4; if f iIf <0, the bow and fire protection system is offline, proceed to step S35.

[0051] S35: Let f i =0,F i =[0 0 0] T And calculate the coordinate vector e of the pantograph-contact network system in the next time step when the system is offline. i+1 The calculation formula is:

[0052]

[0053] Step 4: Determine whether the contact force calculation for all time steps has been completed. If it has, output the pantograph-catenary contact force sequence for all time steps and end the calculation. If it has not been completed, update the pantograph position according to the vehicle speed, update the iterative time step index and pantograph-catenary state, and return to Step 3 to perform the calculation for the next time step.

[0054] The beneficial technical effects of this invention compared to the prior art are as follows:

[0055] This invention employs a method for estimating the relative motion velocity between the pantograph and the catenary, avoiding the iterative steps required to determine the direction of the contact force. The system is linearized in its initial equilibrium state, avoiding the assembly of the system stiffness matrix required for calculating elastic internal forces at each time step. The main matrices are sparsified, and the generalized stiffness matrix is ​​decomposed using LU decomposition, reducing computational complexity. By delaying the contact state determination step, the additional computational load required for contact state determination is avoided. Compared with traditional pantograph-rigid catenary simulation methods, this invention significantly improves system simulation efficiency, substantially reduces simulation time costs, and ensures extremely high computational accuracy. Attached Figure Description

[0056] Figure 1 This is a flowchart of the rapid simulation method for the three-dimensional dynamic contact behavior of the pantograph rigid contact network according to the present invention.

[0057] Figure 2 This is a comparison chart of the results of the fast simulation method proposed in this invention and the traditional method. Detailed Implementation

[0058] The present invention will now be described in further detail with reference to the accompanying drawings and specific embodiments.

[0059] This embodiment presents a simulation case of dynamic contact between a 30-span rigid contact wire and a pantograph on a straight section, with a vehicle speed of 200 km / h. The method proposed in this invention is applicable to various beam elements. In this embodiment, an Absolute Nodal Coordinate (ANCF) beam element model of the rigid contact wire is established. The pantograph-rigid contact force is calculated using both the standard three-dimensional contact force simulation method and the rapid simulation method for three-dimensional dynamic contact behavior of the pantograph-rigid contact wire proposed in this invention. The calculation results of the two methods are then compared. Specific parameters in this embodiment are shown in Tables 1 and 2, with a friction coefficient of 0.3 between the pantograph head and the contact wire. In this embodiment, a right-handed Cartesian coordinate system is used, with the train's forward direction as the positive X direction, the positive Y direction as the vertically upward direction, and the Z direction perpendicular to the track centerline.

[0060] Table 1. Parameters of the rigid contact wire in the embodiments.

[0061]

[0062] Table 2 shows the pantograph parameters in the embodiments.

[0063]

[0064] The process of this invention is as follows Figure 1 As shown, the specific implementation process is as follows:

[0065] Step 1: Based on the rigid contact wire parameters shown in Table 1, establish a finite element model of the rigid contact wire using absolute nodal coordinate beam elements (this invention is applicable to all types of beam elements) and three-dimensional spring elements. Based on the pantograph parameters shown in Table 2, establish a pantograph model using the reduced mass block method, and construct the motion equations of the pantograph-rigid contact wire. The specific implementation method is as follows:

[0066] S11: Treating the busbar and contact wire as a whole, and discretizing it using absolute nodal coordinate beam elements according to the rigid contact network design parameters, the coordinate vectors of each element and node are determined. Based on the absolute nodal coordinate beam element theory, the mass matrix, stiffness matrix, damping matrix, and load vector of each element are calculated.

[0067] S12: The suspension structure and clamps are considered as a whole and simulated using a three-dimensional spring element. The mass matrix M of the three-dimensional spring element is... s and stiffness matrix K s They are respectively:

[0068]

[0069]

[0070] Where, m s k is the equivalent mass of the suspension structure.x k y k z These are the equivalent stiffnesses of the three-dimensional spring on the three coordinate axes X, Y, and Z of the global coordinate system, respectively.

[0071] S13: Based on the parameters shown in Table 2, the mass block model of the pantograph is established using the reduced mass block method. The motion equations of the pantograph three-mass block model are as follows:

[0072]

[0073] Where m1, m2, and m3 are the masses of the bow head, upper frame, and lower frame, respectively; c1, c2, and c3 are the damping of the bow head, upper frame, and lower frame, respectively; k1, k2, and k3 are the stiffness of the bow head, upper frame, and lower frame, respectively; y1, y2, and y3 are the displacements of the bow head, upper frame, and lower frame, respectively; F L This refers to the static lifting force exerted on the pantograph.

[0074] S14: The motion equations of the pantograph-rigid contact network system are constructed using the standard finite element method as follows:

[0075]

[0076] Where M and C are the overall mass matrix and damping matrix of the pantograph-rigid contact network system, respectively; Let Q be the velocity and acceleration vectors of the system as a whole; let Q be the elastic internal force vector of the system as a whole; let P be the external load vector acting on the system; the formula for calculating the elastic internal force vector of the system is:

[0077]

[0078] Where K is the system stiffness matrix and e is the system coordinate vector.

[0079] Step 2: Substitute the gravity load into the external load vector in the pantograph-rigid contact wire motion equation, and use Newton's iteration method to calculate the coordinate vector e0 of the system in the initial equilibrium state under gravity; based on the initial equilibrium state of the pantograph-contact wire system under gravity, pre-calculate the elastic internal force vector Q0 of the system in this state, and assemble the tangential stiffness matrix K of the system in this state. t The generalized stiffness matrix J is then calculated and sparsified. The expression for the system's generalized stiffness matrix J is:

[0080]

[0081] Here, Δt represents the iteration time step in the dynamic simulation process. The zero vector serves as a placeholder; its value needs to be updated in subsequent dynamic simulations. Memory is allocated in advance here to reduce the time overhead of memory reallocation.

[0082] Step 3: Calculate the contact force vector of the pantograph-rigid contact network system at the current time step, and the system coordinate vector at the next time step.

[0083] S31: Calculate the system's elastic internal force Q at the current time step using the concept of linearization. i And calculate the generalized load vector. The calculation formulas are as follows:

[0084] Q i =Q0+K t (e i -e0) (7)

[0085]

[0086] Where the subscript i is the time step index. At the initial time, i.e., when i = 0, e i-1 Calculated by the following formula:

[0087]

[0088] S32: Calculate the contact force direction vector n using the estimated relative velocity. i The calculation formula is:

[0089]

[0090] Where μ is the Coulomb friction coefficient between the bow and the fire net, and n n,i The normal contact force direction vector is...

[0091]

[0092] Among them, A p,i and A c,i These are the tangential vectors of the pantograph head and the contact wire, respectively.

[0093] n t,i The estimated tangential force direction vector.

[0094]

[0095] Among them, v i The estimated relative velocity between the bow head and the contact line at the contact point is given by the following formula:

[0096]

[0097] Among them, v t Let G be the vehicle's velocity vector. Let G be the contact constraint matrix, written as...

[0098]

[0099] Define the pantograph's forward direction as "forward" and the opposite direction as "rear," where n1 is the number of degrees of freedom of all nodes behind the beam elements in contact on the rigid contact network, and n2 is the number of degrees of freedom of all nodes in front of the beam elements in contact on the rigid contact network; S c S is the shape function matrix of a beam element currently in contact. p This is an auxiliary matrix used to calculate the vertical position of the bow head, written as:

[0100]

[0101] S33: Assuming the pantograph-catenary system is in contact, update the last row and last column of the system's generalized stiffness matrix and perform LU decomposition on it. The update format is as follows:

[0102]

[0103] Where Δt is the iteration time step, i is the iteration index, and T = [0 1 0] is the auxiliary vector. Performing LU decomposition on J yields the upper triangular matrix J. U and lower triangular matrix J L Satisfying J L ·J U =J.

[0104] S34: Calculate the contact force at the current time step and the system coordinate vector at the next time step. The calculation formula is as follows:

[0105]

[0106] Among them, U i It includes the system's coordinate vector and the magnitude of the contact force.

[0107]

[0108] f i The total contact force is the magnitude of the total contact force at the current time step. The total contact force vector is calculated by the following formula:

[0109] F i =f i ·n i (19)

[0110] The actual contact state between the pantograph and the rigid contact wire is determined by the sign of the contact force. If f i If f ≥ 0, then the bow and fire net are in contact, proceed to step 4; if f i If <0, the bow and fire protection system is offline, proceed to step S35.

[0111] S35: Let f i =0,F i =[0 0 0]T And calculate the coordinate vector e of the pantograph-contact network system in the next time step when the system is offline. i+1 The calculation formula is:

[0112]

[0113] Step 4: Determine whether the contact force calculation for all time steps has been completed. If it has, output the pantograph-catenary contact force sequence for all time steps and end the calculation. If it has not been completed, update the pantograph position according to the vehicle speed, update the iterative time step index and pantograph-catenary state, and return to Step 3 to perform the calculation for the next time step.

[0114] The contact pressure between the pantograph and the rigid contact wire can be quickly calculated through the above steps. The models in Tables 1 and 2 are simulated using the standard three-dimensional contact force simulation method for pantograph-rigid contact wire, and the results are compared with those obtained by the rapid simulation method proposed in this invention. The total contact force obtained by the two methods is compared as follows: Figure 2 As shown in Table 3, the statistical comparison of the results is presented.

[0115] Table 3 Comparison of results obtained from the two methods

[0116]

[0117] It can be seen that the results obtained by the rapid simulation method are basically consistent with those obtained by the standard calculation method. The absolute value of the relative deviation of the total contact force is only 0.48%, indicating that the rapid simulation method proposed in this invention has almost no impact on the accuracy of the results. The calculation time required by the standard simulation method is 2092.89s, while the calculation time using the rapid simulation method proposed in this invention is 48.75s, which is only 2.33% of the standard method, saving 97.67% of the calculation time cost. The results show that the rapid simulation method for the three-dimensional dynamic contact behavior of pantograph-rigid contact network proposed in this invention can significantly reduce the calculation time cost of simulation without significantly affecting the accuracy of the results.

Claims

1. A method for fast simulation of the three-dimensional dynamic contact behavior of a rigid catenary-pantograph system, characterized in that, Includes the following steps: Step 1: Use beam elements to simulate the busbar and contact wire, use spring elements to simulate the suspension structure and wire clamps, establish a finite element model of the rigid contact network, use the reduced mass block method to establish the pantograph model, and construct the motion equations of the pantograph-rigid contact network. S11: Based on the design parameters of the rigid contact network, the rigid contact network is discretized using beam elements to determine the coordinate vectors of each element and node; based on beam element theory, the mass matrix, stiffness matrix, damping matrix, and load vector of each element are calculated. S12: The suspension structure and the wire clamp are regarded as a whole, simulated by using a three-dimensional spring element, and the mass matrix and the stiffness matrix of the three-dimensional spring element are respectively: and ​ （1） （2） wherein, is the equivalent mass of the suspended structure, , , respectively the equivalent stiffness of the three-dimensional spring on the three coordinate axes of the global coordinate system XYZ. S13: Using the reduced mass block method, a mass block model of the pantograph is established. The motion equations of the pantograph three-mass block model are: （3） in, , , The masses of the bow head, upper frame, and lower frame are respectively. , , These are the damping components for the bow head, upper frame, and lower frame, respectively. , , These are the stiffnesses of the bow head, upper frame, and lower frame, respectively. , , These represent the displacements of the bow head, upper frame, and lower frame, respectively. This refers to the static lifting force exerted on the pantograph. S14: The motion equations of the pantograph-rigid contact network system are constructed using the standard finite element method as follows: （4） in, , These are the mass matrix and damping matrix of the pantograph-rigid contact wire system as a whole; , These are the overall velocity and acceleration vectors of the system, respectively. This represents the overall elastic internal force vector of the system. The vector of external loads acting on the system; the formula for calculating the vector of elastic internal forces in the system is: （5） in, Here is the system stiffness matrix. The system coordinate vector; Step 2: Substitute the gravity load into the external load vector in the pantograph-rigid contact network motion equation, and use Newton's iteration method to calculate the coordinate vector of the system in the initial equilibrium state under gravity; based on the initial equilibrium state of the pantograph-contact network system under gravity, pre-calculate the elastic internal force vector of the system in this state, assemble the system tangential stiffness matrix and generalized stiffness matrix in this state, and perform sparsification on them. Among them, the system's generalized stiffness matrix The expression is: （6） in, This refers to the iteration time step in the dynamic simulation process; Step 3: Calculate the contact force vector of the pantograph-rigid contact network system at the current time step, and the system coordinate vector at the next time step; S31: Calculate the system's elastic internal forces at the current time step using the concept of linearization. And calculate the generalized load vector. The calculation formulas are as follows: （7） （8） Among them, subscript For time step index, Let be the coordinate vector of the system in its initial equilibrium state under the influence of gravity. Let be the elastic internal force vector of the system in this state. The tangent stiffness matrix of the system under this condition; This represents the external load vector at step i; At the initial moment, i.e. hour, Calculated by the following formula: （9） S32: Calculate the contact force direction vector using the estimated relative velocity. The calculation formula is: （10） in, The coefficient of friction between the bow and the fire net is the Coulomb friction coefficient. The normal contact force direction vector is... （11） in, and These are the tangential vectors of the pantograph head and the contact wire, respectively. The estimated tangential force direction vector. （12） in, The estimated relative velocity between the bow head and the contact line at the contact point is given by the following formula: （13） in, The vehicle's speed vector; For the contact constraint matrix, write: （14） Define the direction of pantograph movement as forward and the opposite direction as backward. Let be the number of degrees of freedom of all nodes behind a beam element in contact on a rigid contact network. This represents the number of degrees of freedom of all nodes in front of the beam element in contact on the rigid contact network. Here is the shape function matrix of the beam element that is currently in contact. This is an auxiliary matrix used to calculate the vertical position of the bow head, written as: （15） S33: Assuming the pantograph-catenary system is in contact, update the last row and last column of the system's generalized stiffness matrix and perform LU decomposition on it. The update format is as follows: （16） in, The iteration time step, For iterative index, As an auxiliary vector, for Performing LU decomposition yields the upper triangular matrix. and lower triangular matrix ,satisfy ; S34: Calculate the contact force at the current time step and the system coordinate vector at the next time step. The calculation formula is as follows: （17） in, It includes the system's coordinate vector and the magnitude of the contact force, that is: （18） The total contact force is the magnitude of the total contact force at the current time step. The total contact force vector is calculated by the following formula: （19） The actual contact state between the pantograph and the rigid contact wire is determined by the sign of the contact force. If the bow and fire net are in contact, proceed to step 4; if If the bow and fire net are offline, proceed to step S35; S35: Order , And calculate the coordinate vector of the pantograph-contact network system in the next time step when it is offline. The calculation formula is: （20） Step 4: Determine whether the contact force calculation for all time steps has been completed. If it has, output the pantograph-catenary contact force sequence for all time steps and end the calculation. If it has not been completed, update the pantograph position according to the vehicle speed, update the iterative time step index and pantograph-catenary state, and return to Step 3 to perform the calculation for the next time step.

Images