Boundary-Adapted Numerical Modeling of Flow in a Hydrostatic Leadscrew

: A new method is presented to model and predict the ﬂow ﬁelds of the hydrostatic leadscrews with greater accuracy. It is different from those available methods, in which various bearings are assumed to be equivalent to the screw-nut pair within a pitch by various means. In this new method, a helical coordinate system adapting to the boundaries of the ﬂow ﬁelds is constructed, which makes the screw-nut meshing clearance calculated more accurate. Based on the ﬁnite difference method (FDM), the meshing clearance is discretized into a number of ﬂow ﬁelds, which are created by numerous couples of parallel-plate elements moving relatively along the helicoid. The numerical model is solved in MATLAB, and the analyses about the pressure ﬁelds demonstrate its favorable performances in reﬂecting the actual ﬂow ﬁelds. Furthermore, the simulation results are compared with the experimental values, conﬁrming the feasibility of the proposed method.


Introduction
The hydrostatic leadscrews are widely used for precision machining, and increasingly high demands are placed on them. The numerical studies of the systems can circumvent many obstacles inherent in the experiments, and thus reduce the test costs, which is of great significance for the developments of the prototype machines.
A drawback of the approximation methods of axial projection and tapered bearing substitution is that the influences of torsion and the lead angle are ignored, thus the actual conditions of flow fields are poorly reflected.
With the lead angle and thread angle both covered, the sector approximation method is superior to the others. In the method, there is a one-to-one mapping from the finitedifference elements of the helicoid to the those of the annular sector. Furthermore, there are some approximate rules to be followed between the sums of the elements of the helicoid and those of the sector.
Theoretically, the areas of two bounded surfaces with different shapes may be different even if the sides of them are of equal length, respectively, such as a bounded plane and a bounded curved surface, with the same length on each boundary. In [11], the restrictions were placed only on the lengths of the flow-field boundaries, not the areas, which may make such parameters as the lead angle, thread angle, arc length, torsion, and curvature of each grid element approximated to varying degrees. If such variations occur, then the calculated flow-field pressures would deviate from the actual values. The reason for this is that the pressure values are strongly related to the areas and shapes of the loading surfaces, although the directions of the speed and the pressure at each node remain consistent before and after approximation. As is well known, the denser the meshes, the smaller the deviations caused by the aforesaid approximations, but the according computing costs will increase significantly. Furthermore, there may be a divergence problem of solutions with the increase of the number of grid elements. Lastly, the grid elements may change to adapt to the different discretization methods of the fluid equations and geometric parameters of the system.
Available studies on the flow-field modeling of the hydrostatic and aerostatic leadscrews are quite few. Although the sector approximation method would be superior to the others theoretically, the solutions based on the method were not strictly verified, and the flow fields were not reflected well. In view of the shortcomings of the current methods, a new modeling method, in which a follow-up coordinate system adapting to the boundaries of the flow fields is constructed along the pitch diameter helix of the nut, is proposed in this research, to provide a more convincing prediction of the flow fields.
In this research, the method is implemented and discussed by MATLAB. Compared with the common flow-field simulations in the CFD professional software, the simulations in MATLAB are realized by converting all the information into MATLAB codes, without three-dimensional geometric models of the flow fields. For a right-hand-thread leadscrew, the nut moves downwards along the axial direction when the screw is turned counterclockwise. Point e in Error! Reference source not found.a,b,d is the starting point of the helical coordinate system, point E in Error! Reference source not found.b is any point on the s -axis, point F in Error! Reference source not found.b, e is a point whose distance from point E in the ' y direction is h , and h is the clearance thickness of the gap flow fields. The normal section " nn " of the screw is shown in Error! Reference source not found.c, the section of the gap fields at the initial time is shown in Error! Reference source not found.d, and the relations among velocity components of point F are shown in Error! Reference source not found.e. (c) " nn " of the screw; (d) initial section; (e) velocity relations.
As shown in Figure 1. Modeling diagram based on the helical coordinate system. (a) Schematic diagram; (b) modeling of fluid domain; (c) " nn " of the screw; (d) initial section; (e) velocity relations.b,d, a helical coordinate system " sxy ", composed of the s -axis and a follow-up coordinate system " oxy ", was constructed on the upper helicoid of the nut [19,20]. The s -axis gone from point e along the pitch diameter helix of the helicoid, and " oxy " changed along the s -axis (tangent vector t ) to adapt to the normal and binormal vectors ( n and b ) of the Frenet frame of pitch diameter helix.

Construction of Helical Coordinate System Adapting to the Boundaries of the Flow Fields
A screw-driven hydrostatic leadscrew with a six-turn-thread nut was chosen as the research object to describe the new modeling method, and four oil chambers were assumed to be evenly spaced on each of the four turns of threads in the middle of the nut.
The leadscrew was assumed to be oriented vertically for the convenience of description. The flow fields on the upper side of the nut were named the upper flow fields, and those on the other side were named the lower flow fields. Besides, the research was carried out under the following assumptions: (1) the pressure gradient along the film thickness direction is 0; (2) the outer boundary pressure is 0; (3) the velocity gradient in the helical direction can be ignored, compared with the thickness direction; (4) the volume force and inertia force can be neglected, compared with the viscous force; (5) the research involves only laminar motion of an incompressible Newtonian fluid; (6) the boundary velocities satisfy the no-slip condition; and (7) the normal section of fit clearance is assumed as rectangular, with a trapezoidal thread discussed in this research.
The Lagrange method was adopted to describe the flow fields in this research [18]. The upper flow fields were taken for example. As shown in Figure 1a,b, "Oξ 1 ξ 2 ξ 3 " is a global coordinate system ("Orθξ 3 " is the corresponding cylindrical coordinate system), and "O 1 XYZ" is a local inertial frame fixed on the nut. The O 1 Z-axis coincided with the Oξ 3 -axis, and the O 1 X-axis intersected the upper helicoid of the nut at point e on the pitch diameter helix. With counterclockwise rotations defined as positive, the screw rotation angle θ was calculated, starting from the position of the O 1 X-axis at the initial time. For a right-hand-thread leadscrew, the nut moves downwards along the axial direction when the screw is turned counterclockwise. Point e in Figure 1a,b,d is the starting point of the helical coordinate system, point E in Figure 1b is any point on the s-axis, point F in Figure 1b, e is a point whose distance from point E in the y direction is h, and h is the clearance thickness of the gap flow fields. The normal section "n-n" of the screw is shown in Figure 1c, the section of the gap fields at the initial time is shown in Figure 1d, and the relations among velocity components of point F are shown in Figure 1e.
As shown in Figure 1b,d, a helical coordinate system "sy", composed of the s-axis and a follow-up coordinate system "oxy", was constructed on the upper helicoid of the nut [19,20]. The s-axis gone from point e along the pitch diameter helix of the helicoid, and "oxy" changed along the s-axis (tangent vector t) to adapt to the normal and binormal vectors (n and b) of the Frenet frame of pitch diameter helix.
The origin o of the follow-up coordinate system could be any point on the pitch diameter helix. R f was employed to locate point o in "O 1 XYZ", and the position of the point in "Oξ 1 ξ 2 ξ 3 " could be expressed as Equation (1): The position of point O 1 in "Oξ 1 ξ 2 ξ 3 " was described by R O 1 , whose relationship with the axial displacement of the nut Z nut was: i, j and k are the unit vectors in the ξ 1 , ξ 2 and ξ 3 direction, respectively. An arbitrary point in the fields could be located in "Oξ 1 ξ 2 ξ 3 " through Equation (3): R hA is the position of the arbitrary point in "sxy" or "O 1 XYZ". "Oξ 1 ξ 2 ξ 3 " and "O 1 XYZ" were coincident at the coordinate origin at the initial moment, that was to say, The Frenet frame of the s-axis in "O 1 XZ" could be solved based on Equations (4)- (6), to evaluate the orthogonality of the helical coordinate system.
The arc length was calculated: Here, P, κ, τ, and r m are the pitch, curvature, torsion, and pitch radius of the helix along the s-axis, respectively. r m was supposed to be equal to the pitch radius of the gap fields. The helical coordinate system "sxy" was considered to be non-orthogonal, according to the judgment rule [21].
Shown in Figure 1b,d, according to the methods for rotating the Frenet frame in [21,22], each follow-up coordinate system "oxy" was rotated by ϕ (ϕ is half of the thread angle of the normal section at the position with a radius of r m in Figure 1. Modeling diagram based on the helical coordinate system. (a) Schematic diagram; (b) modeling of fluid domain; (c) "-" of the screw; (d) initial section; (e) velocity relations.c) around the s-axis to obtain a new follow-up coordinate system "ox y ". In the corresponding new helical coordinate system "sx y ", the normal and binormal vectors (n and b ) of the new Frenet frame were aligned with the corresponding flow field boundaries, respectively.
The normal and binormal vectors before and after rotation satisfied the relationships: The position of any selected point in "sx y " after rotation could be described as: Similar to Equation (7), the arc length was calculated: It could be deduced from Equation (11) that the new helical coordinate system was not orthogonal until τ = ∂ϕ ∂s , and then there was such a relation as follows: In light of the above analyses, the helical coordinate system was deduced to be nonorthogonal with a constant ϕ.
As mentioned above, the Lagrangian description was employed during the process. The relative motion of the fluid particle in "O 1 XYZ" was analyzed firstly for convenience, and then the absolute motion in "Oξ 1 ξ 2 ξ 3 " was solved based on the theory of relative motion (refer to the relations shown in Figure 1e).

Derivation of Reynolds Equation in Helical Coordinate System
In the following analyses, the group of symbols "s, x, y" was utilized to replace "s, x , y " after rotation, and to locate the flow-field particles in the new helical coordinate system, for ease of description. As mentioned in Section 2.1, the motions of fluid particles relative to the nut are analyzed in this section.
According to the assumptions in (1)- (7) in Section 2.1, the Navier-Stokes equation could be simplified as follows [21][22][23]: The continuity equation could be reduced as follows: where M = 1 − κx cos ϕ − κy sin ϕ, and the velocities of any fluid particle in the x, y, and s direction were indicated by such vectors as "u, v, w", respectively. M = 1 − κx cos ϕ was assumed, as the order of magnitude in the y direction was much smaller than that in other directions. Then, the first formula in Equation (14) was developed into: It could be deduced that ∂p ∂y = 0 from the assumption in (1) in Section 2.1. The upper flow fields were taken for example. The fluid particle, which was stationary relative to the nut, was assumed to have a y-coordinate of 0 ("u 1 , v 1 , w 1 "). There was no slip motion between the screw and the particle whose y-coordinate was h ("u 2 , v 2 , w 2 ").
When Equation (16) was integrated twice with respect to y, ω was obtained: where, for x ∈ ∀, w(x) is true as long as y = 0. When the second formula in Equation (14) was integrated twice with respect to x, u was given by: The unit discharges in the s and x direction could be solved by: Equation (15) was developed into: Here, the term y ∂ω ∂x could be ignored, because ∂ω ∂x << ∂ω ∂y . When Equation (21) was integrated with respect to y, the simplified Reynolds equation [24,25] was obtained:

Boundary Conditions of the Helical Flow Fields
In the light of the kinematics principle of particle and the theory of moving coordinate system [26,27], the relations among the velocity components of the fluid particle, which makes no slip motion relative to the screw, are shown in Figure 1e. Then, a set of relations could be deduced from v Here, v r → v r is the relative velocity of a fluid particle, v x , v y , and v s are the velocities of the particle in the x, y, and s direction, respectively, δ is the angle between v r and v s , and • Z nut is the axial velocity of the nut. And then, where λ r is the lead angle of the helix at the position with a radius of r, ∆h is the transient variation of clearance thickness, ∆Z is the transient displacement from the equilibrium position of the nut caused by ∆h, and ∆ • Z is the change rate of ∆Z. Thus, v s could be defined by: When v x was ignored, v y ≈ v r sin δ. Then, the boundary conditions of the flow fields are presented in Table 1

Lower Fields
Parameter Value Parameter Value

Discretization of Reynolds Equation
The following relations were supposed: Here, ∆h s is the steady component of ∆h, [28]. Through the definitions: , the dimensionless transient Reynolds equation of the upper flow fields was given by: The five-point finite difference method (FDM) was employed to discretize the Reynolds equation in Equation (33), and the kth iterative pressure was given by: The steady state were taken as an example. The coefficients are given by: where "i" and "j" are applicable to all internal nodes.

Solutions of Pressure and Flow
As shown in Figure 2, fluids were assumed to flow out of an oil chamber in four directions, and the distribution of oil chambers in one of the four turns of helical fields was obtained by means of an axial projection [29]. The chambers were counterclockwise named as "Chamber 1", "Chamber 2", "Chamber 3", and "Chamber 4".
Here, Q out is the total flow of an oil chamber for a certain time.
Here, out Q is the total flow of an oil chamber for a certain time.
Continuing from Section 2.4.1, the upper flow fields in the steady state were still taken for example. On the basis of Equation Error! Reference source not found. and Error! Reference source not found., the dimensionless discharges in each direction could be given as follows: 3 Figure 2. Axial projection of the chamber distributions.
Continuing from Section 2.4.1, the upper flow fields in the steady state were still taken for example. On the basis of Equations (19) and (20), the dimensionless discharges in each direction could be given as follows: Q out could be solved by summing over Equations (42)-(45).
The capillary restrictor was taken as an example. The relationship between the chamber pressure and flow could be obtained according to the flow continuity through an oil chamber: here, Csc is the throttle coefficient. And the node pressures in the same oil chamber were considered to be equal. Based on Equations (34)-(46), the node pressures of the upper fields in the steady state could be solved through the method of successive over-relaxation (SOR) [28,[30][31][32].

Solution of Static and Dynamic Characteristic Parameters
Only the axial static and dynamic characteristics were studied in this research.
where W is the axial bearing capacity, S is the axial stiffness coefficient, and C is the axial damping coefficient. While W > 0, the resultant force applied on the nut was upward, otherwise, it was downward. Additionally, h Td and h Bd are the film thicknesses of the upper and lower flow fields in the disturbed state, respectively. As shown in Equations (47)-(49), the effects of viscous force were considered in this research.

Flow-Field Pressure Analysis
The main geometric parameters of the system are listed in Table 2, and the numerical model of the flow fields is solved by MATLAB software. The negative pressure values outside film-rupture boundary which is taken as 0 MPa [33], are set to zero, when the steady pressures are solved with MATLAB, but it is not necessary for those under the disturbed state to be processed in that way. It can be found in the calculation that the first-order pressure disturbances are solved based on the steady pressures, and the total pressures are less influenced by the second-order pressure disturbances than the others. Thus, the steady pressures can be taken as the main targets, when the pressure fields are analyzed.
Dimensionless pressure distributions of the flow fields based on the new method, are shown in Figure 3. The number of grid points in the x direction is 100, and that in the s direction is 360. Ω is 500 rm, is 0.4, "X, Y, Z" is used to orientate the fluid particle, "T" represents the upper fields of the nut, and "B" represents the lower fields, with "Z d " and " • Z d " indicating the first-order and second-order pressure disturbance, respectively, for example, "P Ts " denotes the dimensionless pressure of the upper flow fields of the nut in the steady state.
As mentioned in Section 2.1, the system is assumed to be oriented vertically, and the sealing surfaces of the upper and lower fields are named the upper sealing surfaces and the lower sealing surface, respectively. Besides, for the flow fields on either side of the nut, the sealing surfaces at the both ends of the four turns of helical fields with oil chambers are called the front sealing surface and the rear sealing surface, respectively, according to their locations relative to the coordinate origin O 1 .
Actuators 2021, 10, x FOR PEER REVIEW 13 of 20 As mentioned in Section 2.1, the system is assumed to be oriented vertically, and the sealing surfaces of the upper and lower fields are named the upper sealing surfaces and the lower sealing surface, respectively. Besides, for the flow fields on either side of the nut, the sealing surfaces at the both ends of the four turns of helical fields with oil chambers are called the front sealing surface and the rear sealing surface, respectively, according to their locations relative to the coordinate origin O1. It can be concluded from Error! Reference source not found.a,d that the nut is subjected to a downward driving force from the fluids, which is consistent with the actual condition mentioned in Section 2.1. Such a result is related to the regulating functions of the restrictors. It can be concluded from Figure 3a,d that the nut is subjected to a downward driving force from the fluids, which is consistent with the actual condition mentioned in Section 2.1. Such a result is related to the regulating functions of the restrictors.
As shown in Figure 3a, a symmetrical drop around r ci ≤ r ≤ r co occurs in the four turns of helical fields with oil chambers, and the pressure values in the front and rear sealing surfaces are close to zero. The flow-field pressures in Figure 3d show irregular change trends with variable gradients in the whole fields, especially in the front and rear sealing surfaces. Such irregular changes also obviously occur in the front and rear sealing surfaces in Figure 3c,e,f. The seemingly irregular changes mainly serve as a result of a combined effect of the torsion, lead angle, thread angle of the helix, and the centrifugal forces, and the effect is larger than that of the restrictors, especially when the fluids are squeezed (∆Z = 0.4). As far as the upper flow fields are concerned, the enlarged clearance results in weaker diffusions to the sealing surfaces, hence the pressures in the sealing surfaces in Figure 3a are too small to show the pressure gradients caused by the above factors.
It is found that the illustrations in Figures 4 and 5 are not intuitive enough, and some less obvious changes are not easy to catch. Thus, as shown in, the variations of the steady pressures with ∆Z in the x(r) and θ(s) direction are solved to observe and analyze the pressure fields more intuitively.
Actuators 2021, 10, x FOR PEER REVIEW 14 of 20 As shown in Error! Reference source not found.a, a symmetrical drop around ≤ ≤ ci co r r r occurs in the four turns of helical fields with oil chambers, and the pressure values in the front and rear sealing surfaces are close to zero. The flow-field pressures in Error! Reference source not found.d show irregular change trends with variable gradients in the whole fields, especially in the front and rear sealing surfaces. Such irregular changes also obviously occur in the front and rear sealing surfaces in Error! Reference source not found.c,e,f. The seemingly irregular changes mainly serve as a result of a combined effect of the torsion, lead angle, thread angle of the helix, and the centrifugal forces, and the effect is larger than that of the restrictors, especially when the fluids are squeezed As far as the upper flow fields are concerned, the enlarged clearance results in weaker diffusions to the sealing surfaces, hence the pressures in the sealing surfaces in Error! Reference source not found.a are too small to show the pressure gradients caused by the above factors.
It is found that the illustrations in Error! Reference source not found. are not intuitive enough, and some less obvious changes are not easy to catch. Thus, as shown in Error! Reference source not found. and Error! Reference source not found., the variations of the steady pressures with ΔZ in the x ( r ) and θ ( s ) direction are solved to observe and analyze the pressure fields more intuitively.  There is a mixture of Couette flow and Poiseuille flow in this research [34]. The motion of the screw relative to the nut brings the Couette flow. Poiseuille flow exists not only between a chamber and the adjacent sealing surface, but also between adjacent oil chambers. When the clearance thickness is small, the Couette flow is dominant to maintain the pressures constant among different oil chambers, such as when in Error! Reference source not found.b. Then, these chamber pressures do not vary in the helical direction until the clearance thickness increases to large enough. As is shown in Error! Reference source not found.a, , not only are pressure differences between the oil chambers and the sealing surfaces shown, but also the pressure differences caused by such factors as torsion and the lead angle are reflected among different oil chambers, which demonstrates the effect of Poiseuille flow on the chamber pressures.
, the reason why the pressures in the sealing surfaces in Error! Reference source not found.a are so small has been explained in Error! Reference source not found.. The fluids on the lower side of the nut are squeezed and then diffuse to the sealing surfaces, which is intensified with the increase of ΔZ . As shown in the sealing surfaces in Error! Reference source not found.b, a few sharp increases in the yellow curve occur when ΔZ is large enough. There is a mixture of Couette flow and Poiseuille flow in this research [34]. The motion of the screw relative to the nut brings the Couette flow. Poiseuille flow exists not only between a chamber and the adjacent sealing surface, but also between adjacent oil chambers. When the clearance thickness is small, the Couette flow is dominant to maintain the pressures constant among different oil chambers, such as when ∆Z = 0.1, 0.3 in Figure 4a and ∆Z = 0.1, 0.3, 0.5 in Figure 4b. Then, these chamber pressures do not vary in the helical direction until the clearance thickness increases to large enough. As is shown in Figure 4a, when ∆Z = 0.5, not only are pressure differences between the oil chambers and the sealing surfaces shown, but also the pressure differences caused by such factors as torsion and the lead angle are reflected among different oil chambers, which demonstrates the effect of Poiseuille flow on the chamber pressures.
When ∆Z = 0, i.e., ∆h = 0, the reason why the pressures in the sealing surfaces in Figure 4a are so small has been explained in Figure 3. The fluids on the lower side of the nut are squeezed and then diffuse to the sealing surfaces, which is intensified with the increase of ∆Z. As shown in the sealing surfaces in Figure 4b, a few sharp increases in the yellow curve occur when ∆Z is large enough.
The chambers contribute to the bearing capacity much greater than the sealing surfaces, which is achieved with the help of the restrictors. As is shown in Figure 4a,b, with the increase of ∆Z, the chamber pressures of the upper flow fields gradually increase, yet the lower ones decrease, which brings about a larger bearing capacity.  Figure 5b,d,f, the positions of fluid particles are transformed from "sxy" to "Orθξ 3 " for more intuitive analyses.
According to Figure 2, "j" denotes a set of nodes in the s direction. "hssB" is the midpoint of pitch diameter helix of the front sealing surface, " q1a" is the intersection a of pitch diameter helix of the front sealing surface and the adjacent "Chamber 1", "q1b" is the node b of "Chamber 1" in Figure 2 "q1a − b" is the node a − b of "Chamber 1" in Figure 2, "q4b" is the intersection b of pitch diameter helix of the rear sealing surface and the adjacent "Chamber 4", and "hssT"is the midpoint of pitch diameter helix of the rear sealing surface. The nodes at "j" denote all the particles which have a s-coordinate of "j" in the following article.
Compared to the others, the pressure values at "q1a" are close to those at "q4b". For example, when ∆Z = 0.1, the order of magnitude of the upper pressure differences between "q1a" and "q4b" reaches 10 −7 , leaving only a curve visualized in Figure 5a,b. Furthermore, the lower pressure values at "q1a", "q4b", "hssB", and "hssT" are close to each other. For the same ∆Z, only one of the curves at "hssB"and "hssT" is visible in Figure 5e,f. Consequently, when the lower chamber pressures are analyzed, the pressures at "q1a − b" and "q1b" can be selected as the main research objects.
As shown in Figure 5b, the pressures remain relatively stable in the chambers, but descend to the sealing surfaces in the ranges of r i ≤ r ≤ r ci and r co ≤ r ≤ r o . When ∆Z = 0.1, 0.3, there are similar situations at "q1a − b" and "q1b" in Figure 5d. With the growing of ∆Z, the descending slope of the pressures in the range of r co ≤ r ≤ r o in Figure 5d becomes smaller. And when ∆Z is large enough, a pressure jump occurs. For example, there is an obvious jump in the range of r co ≤ r ≤ r o at "q1b" in Figure 5d, when ∆Z = 0.5.
As mentioned in Figures 3 and 4, the bearing capacity of the system mainly benefits from the regulations of the restrictors on the upper and lower chamber pressures. Furthermore, the pressures are also influenced by the torsion, lead angle and thread angle of the helix, and the centrifugal forces. Under the effect of these factors, the diffusions to the sealing surfaces are enhanced when the fluids in the gap are squeezed. As shown in Figure 5c,d, the regulation effects of the restrictors on the lower pressures are weakened by the diffusions. Meanwhile, with the enlarged clearance and the diffusions to the sealing surfaces, larger upper pressures are caused by the restrictors in Figure 5a,b. With the increase of ∆Z, the regulation effects of the restrictors are enhanced, then smaller decline ranges of the lower pressures in Figure 5c,d serve as a result of the combined effect of the restrictors, the torsion, lead angle and thread angle of the helix, and the centrifugal forces. In this situation, the increasing amplitudes of the upper pressures become larger in Figure 5a,b, which reflects the sensitivity of bearing capacity to ∆h. Such a phenomenon is consistent with that mentioned inFigure 4.
The pressure values in the upper sealing surfaces are near zero according to Figure 3, thus they are not discussed in this part. As mentioned above, the pressures in the lower sealing surfaces increase when the fluids are squeezed; such changes are shown in Figure 3d-f, Figure 5e,f, under the combined influence of the torsion, lead angle and thread angle of the helix, and the centrifugal forces. The increasing amplitudes of the pressures become larger with the growing of ∆Z, as shown in Figure 5e,f. The reason why the pressures decrease to zero in the ranges of r i ≤ r ≤ r ci and r co ≤ r ≤ r o is that the fluids diffuse to the root and crest, and the tops of the hills in Figure 5f are found to be closer to the position with a radius of r co , largely due to the centrifugal forces. Such explanations can also be applicable to the conditions at "q1a" and "q4b" in Figure 5d.

Feasibility Evaluation Based on the Bearing Capacity
The aerostatic leadscrew system in [6] is solved numerically based on the method proposed in this research, and the simulation results are compared with the experimental data [6] to evaluate the feasibility of the new method.
The comparisons between the simulation results and experimental data are shown in Figure 6, in which the experimental data for a two-pitch travel are represented by the yellow , and the simulation results at Ω = 120 rpm, Ω = 100 rpm, and Ω = 60 rpm are represented by the red , blue •, and green , respectively. The aerostatic leadscrew system in [6] is solved numerically based on the method proposed in this research, and the simulation results are compared with the experimental data [6] to evaluate the feasibility of the new method.
The comparisons between the simulation results and experimental data are shown in Error! Reference source not found., in which the experimental data for a two-pitch travel are represented by the yellow ★, and the simulation results at = 120rpm Ω , = 100rpm Ω , and = 60rpm Ω are represented by the red ▲, blue •, and green ◆, respectively. It is found that the numerical predictions agree well with the experimental data qualitatively. However, from a quantitative point of view, there are still some deviations between them. It is inevitable that the data points in the two cases would not coincide. Because the information about the speeds were not provided, the experimental results were not collected until the system stopped running and kept steady, rather than real-time like the simulation results. Besides, it is also necessary to consider the divergences caused by the experimental environments, operating methods, geometric errors, etc.
As shown in Error! Reference source not found., when , the larger the value of Ω , the closer the simulation results are to the experimental ones, however, the situation is the opposite in the range of < Δ 1.4μm Z , which is closely associated with the collection method of the experimental data. Among the simulation data, the red ▲ data are closer to the yellow ★ ones with a maximum difference of about 1.8N , which is approximately 14% of the corresponding experimental value.
Therefore, the new method can be considered to be feasible, and the results can provide references for predicting the actual characteristics of the system.

Conclusions
In this research, a new method is attempted to overcome the shortcomings of the existing modeling methods of the hydrostatic and aerostatic leadscrews, and a screwdriven hydrostatic leadscrew is selected as the main research object to describe the proposed method. It is found that the numerical predictions agree well with the experimental data qualitatively. However, from a quantitative point of view, there are still some deviations between them. It is inevitable that the data points in the two cases would not coincide. Because the information about the speeds were not provided, the experimental results were not collected until the system stopped running and kept steady, rather than real-time like the simulation results. Besides, it is also necessary to consider the divergences caused by the experimental environments, operating methods, geometric errors, etc.
As shown in Figure 6, when 1.4 µm ≤ ∆Z ≤ 4.8 µm, the larger the value of Ω, the closer the simulation results are to the experimental ones, however, the situation is the opposite in the range of ∆Z < 1.4 µm, which is closely associated with the collection method of the experimental data. Among the simulation data, the red data are closer to the yellow ones with a maximum difference of about 1.8N, which is approximately 14% of the corresponding experimental value.
Therefore, the new method can be considered to be feasible, and the results can provide references for predicting the actual characteristics of the system.

Conclusions
In this research, a new method is attempted to overcome the shortcomings of the existing modeling methods of the hydrostatic and aerostatic leadscrews, and a screw-driven hydrostatic leadscrew is selected as the main research object to describe the proposed method.
In the new method, the follow-up coordinate system is adjusted to adapt to the boundaries of the flow fields, which makes the calculations of h and ∆Z more accurate. Furthermore, compared with those approximations in the cylindrical coordinate system, the description in the helical coordinate system can provide a better reappearance of the actual shape of the flow fields.
Through the analyses about the pressure fields in Section 3.1, it can be found that the influences of the torsion, lead angle, and thread angle of the helix are taken into consideration in this research, and the actual flow fields are reflected well by the simulation results based on the new method.
The experiments on an aerostatic leadscrew in [6] are employed to discuss the feasibility of the new method in Section 3.2. In the numerical calculations, known parameters in [6] are set to keep the consistency between the conditions of the experiment and the simulations, and Ω is set in different ranges to reduce the interferences caused by the data collection method. In Figure 6, the red curve shows a better consistency with the yellow target curve than the others, and the relative errors between them are mostly within 10%, demonstrating that the new modeling method is possessed of good practicability and generality.
The method in this research is suitable for all helical gap flow fields created by a screwdriven screw-nut pair, and some adjustments may occasionally be made on the assumptions in (1)- (7) or the fluid governing equations in Section 2.1 to satisfy different lubrication media or applications. Although the assumptions are made for the sake of simplicity in this research, the proposed method provides a new way of thinking for the flow-field modeling of various leadscrews with fluids as lubrication media. With the adjustments of the flow-field models, not only various types of fluid-lubricated leadscrews, but also more complicated flow-field properties and phenomena can be solved and analyzed in future researches.

Conflicts of Interest:
The authors declare no conflict of interest.

Nomenclature
Oξ 1 ξ 2 ξ 3 global coordinate system Orθξ 3 cylindrical coordinate system corresponding to Oξ 1 ξ 2 ξ 3 O 1 XYZ local coordinate system on the nut oxy/ox y follow-up coordinate system sxy/sx y helical coordinate system i, j, k unit vectors in the ξ 1 , ξ 2 and ξ 3 direction r r-coordinate in Orθξ 3 r i inner radius of the nut r o outer radius of the screw r m pitch radius of the gap fields r ci inner radius of the chamber r co outer radius of the chamber κ curvature τ torsion ε dimensionless curvature g dimensionless torsion