Development of A Linear Delta Robot with Three Horizontal-Axial Pneumatic Actuators for 3-DOF Trajectory Tracking

: This paper focuses on developing a pneumatic-driven and horizontal-structure linear delta robot (PH-LDR) and increasing its trajectory tracking performance on the three-degrees-of-freedom (3-DOF) space. With the investigation of inverse and forward kinematics, the parallel mechanism of PH-LDR is designed by using three high-power and low-cost rod-less pneumatic actuators (PAs) to track 3-DOF motion, and this has a horizontal structure to enlarge the workspace. Since the PH-LDR features nonlinear coupling among its three axes and is disturbed by the three high-nonlinear rod-less PAs, the tracking control performance is signiﬁcantly decreased, subject to uncertain nonlinearity and parametric uncertainty. Therefore, a fuzzy-PID controller is used to achieve highly accurate 3-DOF trajectory tracking, and furthermore, this study exploits neural networks (NNs) to pre-compensate the impacts arising from the compressibility of air and temperature change. The control system for the PH-LDR also features an embedded controller that allows real-time control. Experimental demonstration veriﬁes the developed PH-LDR with the proposed controller, and the dynamic tracking accuracy in 3-DOF trajectory can be achieved.


Introduction
Linear delta robots (LDRs) for industrial applications are superior to series manipulators, because they have a large carrying capacity and are highly precise and fast [1,2], so they are widely used for assembly detection, pick-and-place, positioning, and packaging [3]. They are constructed using multiple independent closed-loop mechanisms, each of which has an arm that connects a fixed base to a moving platform. LDRs require three linear actuators to drive the arms and to translate one-dimensional linear motion into three-degrees-of-freedom (3-DOF) movements. Linear actuators can be installed on or near to a fixed base, so that most LDRs use lightweight material for the arms and have a high rigidity-to-weight ratio. LDRs also feature excellent positioning accuracy, because the position error from one single closed-loop kinematic chain is spread among the other chains. However, LDRs have some intrinsic drawbacks, such as a small workspace and complicated kinematics, but there has been increased demand for them in recent years.
Pneumatic actuators (PAs) are low-cost and widely used devices for manufacturing. They use compressed air as the operating fluid to push a piston and produce linear movement for payloads [4][5][6]. PAs are highly suited to industrial environments, because they are clean, safe, and easier to maintain and have a high power-to-weight ratio. However, PAs are less accurate in terms of position control and less versatile and flexible for engineering applications than electrical motors with the same power, due to the compressibility of air, air leakage, a low natural frequency and high nonlinearity. Despite these inherent drawbacks, PA-driven facilities are still commonly used in manufacturing.
Milutinovic et al. [7] presented three different LDRs with actuated arms at specific angles to each other: a vertical type, a tri-pyramid type and a horizontal type. For a specific sensitivity, proper design for LDRs can fulfil a larger workspace, but cause more difficult control problems in real time, because of the multi-input and multi-output system dynamics involving various nonlinearities and coupling effects. Choi et al. [8] presented a hybrid serial-parallel pneumatic-driven robot with a high ratio of strength-to-moving-weight for mounting window glass and fixing pixels. Chiang et al. [9] developed a pneumatic-driven LDR with three pneumatic actuators on three vertical axes at 120 degrees to each other, to achieve complex 3-DOF motion control. Furthermore, Chiang et al. [10] developed a three-axial tripod-type pneumatic parallel LDR. A geometric approach was used to derive the kinematic model and drive the LDR in a 3-DOF space with an acceptable output tracking error. In [11], Andrzej et al. designed a 3-DOF tripod-type pneumatic LDR for municipal waste recycling that demonstrated excellent accuracy and repeatability in pick-and-place tasks. Lafmejani et al. [12] presented a 6-DOF pneumatic Gough-Stewart parallel robot. In order to determine an accuracy model for the 6-DOF pneumatic robot, a genetic algorithm (GA) was used to identify unknown parameters within the dynamic model of the cylinder and the proportional valve. Yan et al. [13] presented a multi-objective optimization methodology for the LDRs that optimizes dexterity, stiffness and space utilization, and showed that the methodology gives optimal solutions that are better than the results using current optimization methodology.
To address the problems of high nonlinearity, low accuracy and low robustness in LDRs, numerous control strategies have been proposed. Dongtao [14] determined the kinematic reliability and sensitivity of LDRs. Fu et al. [15] designed an LDR with three arms and showed that problems with kinematic accuracy are mostly due to U joint errors, clearances and driving errors. When the dynamics of a PA are linearized at a fixed operating point to determine the linear state equation, Liu and Bobrow [16] used a linear optimal controller to regulate the pneumatic servo system for testrobotic applications. Richer et al. [5] and Bone et al. [17] proposed a sliding mode control (SMC) to control the motion for Pas, because the SMC compensates for nonlinearities. A SMC uses a model-based approach, which can produce an unstable system which is unstable due to modeling error. To attenuate the modeling error effect in the tracking control performance, Chiang et al. [9,10,18] proposed an adaptive SMC that regulates the position of the moving platform of LDR. Lee et al. [19,20] employed a function approximation technique to determine an unknown model for PAs and increase tracking performance over time for the overall system. Despite the pervasive development of control algorithms in literature and the widespread using of PAs, the practical solution to achieve positioning accuracy and repeatability still remains an open problem in pneumatic-driven systems.
This study presents the development and implementation of a pneumatic-driven and horizontal-structure linear delta robot (PH-LDR) for industrial manufacturing that is cost-effective and highly reliable. The PH-LDR is a pneumatic-driven and horizontal-structure (PH) linear delta robot (LDR), that uses pneumatic actuators to achieve high trajectory tracking performance and reduce the cost, and which is designed as a horizontal structure to enlarge the workspace. The PH-LDR uses three rod-less PAs as the linear actuator to produce the 3-DOF motion, as well as adopting a hybrid control technique to compensate for nonlinearities and uncertainties. The rod-less PAs possess a long stroke, because unlike general pneumatic actuators, the stroke is contained within the envelope of a cylinder body. Without the need of parameter identification of the pneumatic actuating system for the controller implementation, the hybrid controller uses a simple and straightforward approach and contains two parts: a fuzzy-PID controller, as designed by using expert knowledge, to achieve the positional accuracy, and a neural network pre-compensator (NNPC), that predicts variations caused by the environmental changes such as temperature, and works with the fuzzy-PID controller. Experiments with the implemented full-scale test rig are conducted to confirm the validness of the PH-LDR for complex 3-DOF trajectory tracking problems.
The remainder of this paper is organized as follows. In Section 2, the layout of the test rig and its experimental setup are described. Furthermore, a kinematics analysis of the PH-LDR is analyzed in Section 3. The controller design for the PH-LDR is given in Section 4. Two experiments of the PH-LDR are presented in Section 5 to investigate the trajectory tracking performance. Finally, Section 6 concludes this paper. Figure 1 illustrates the test-rig layout for the PH-LDR. The PH-LDR uses three rod-less PAs that are installed under the three links at an angle of 120 to each other. Each of the rod-less PAs produces a linear movement that drives the arms, so that the moving platform moves with 3-DOF. The rod-less PAs are Festo DGC-25-500-KF-PPV-A with a piston that has a 500 mm long stroke and a 25 mm diameter. A proportional servo valve (Festo MPYE-5-1/8-HF-010-B) controls the compressed air that flows into two chambers. The air acts on the surface of the piston to produce a rectilinear translation. The piston remains at the center of the rod-less Pas, while both of the two chambers have the same air pressure of six bars if the proportional servo valve is set to 5 V, as shown in Figure 2. Figure 1 shows that two pressure sensors (Festo SDE1) are installed on the air inlets of the cylinders to measure the pressure inside the two chambers, and three optical linear encoders with a resolution of 1µm are installed on the links to measure the position of the pistons. An embedded controller c-RIO is embedded in the control loop to achieve real-time control. It receives all measurement data via interface cards and sends an input voltage command to proportional servo valves. The interface cards include a NI-9411 digital I/O converter, a NI-9263 analog output interface card, and a NI-9215 analog input interface card. An industrial PC provides a LabVIEW integrated development environment (IDE) for developers and downloads programs into the c-RIO.

Test Rig Layout
Appl. Sci. 2020, 10, x FOR PEER REVIEW  3 of 21 controller. Experiments with the implemented full-scale test rig are conducted to confirm the validness of the PH-LDR for complex 3-DOF trajectory tracking problems. The remainder of this paper is organized as follows. In Section 2, the layout of the test rig and its experimental setup are described. Furthermore, a kinematics analysis of the PH-LDR is analyzed in Section 3. The controller design for the PH-LDR is given in Section 4. Two experiments of the PH-LDR are presented in Section 5 to investigate the trajectory tracking performance. Finally, Section 6 concludes this paper. Figure 1 illustrates the test-rig layout for the PH-LDR. The PH-LDR uses three rod-less PAs that are installed under the three links at an angle of 120 to each other. Each of the rod-less PAs produces a linear movement that drives the arms, so that the moving platform moves with 3-DOF. The rodless PAs are Festo DGC-25-500-KF-PPV-A with a piston that has a 500 mm long stroke and a 25 mm diameter. A proportional servo valve (Festo MPYE-5-1/8-HF-010-B) controls the compressed air that flows into two chambers. The air acts on the surface of the piston to produce a rectilinear translation. The piston remains at the center of the rod-less Pas, while both of the two chambers have the same air pressure of six bars if the proportional servo valve is set to 5 V, as shown in Figure 2. Figure 1 shows that two pressure sensors (Festo SDE1) are installed on the air inlets of the cylinders to measure the pressure inside the two chambers, and three optical linear encoders with a resolution of 1μm are installed on the links to measure the position of the pistons. An embedded controller c-RIO is embedded in the control loop to achieve real-time control. It receives all measurement data via interface cards and sends an input voltage command to proportional servo valves. The interface cards include a NI-9411 digital I/O converter, a NI-9263 analog output interface card, and a NI-9215 analog input interface card. An industrial PC provides a LabVIEW integrated development environment (IDE) for developers and downloads programs into the c-RIO.

Analysis of Kinematics
This paper uses a vector-loop closure equation [21] for each arm to solve the solutions of the inverse kinematic model for the PH-LDR. Figure 3 shows the coordinate frames for the PH-LDR. The three arms of the PH-LDR have an identical kinematic structure. Two coordinate frames are used to analyze the kinematic model of the PH-LDR. There is a fixed Cartesian reference frame { O }, for which the origin is located at the point where T 0 0 cos(( 1) 120 ) sin(( 1) 120 ) 0

Analysis of Kinematics
This paper uses a vector-loop closure equation [21] for each arm to solve the solutions of the inverse kinematic model for the PH-LDR. Figure 3 shows the coordinate frames for the PH-LDR. The three arms of the PH-LDR have an identical kinematic structure. Two coordinate frames are used to analyze the kinematic model of the PH-LDR. There is a fixed Cartesian reference frame {O}, for which the origin is located at the point O(x 0 , y 0 , z 0 ) on the fixed base. The frame {O} also serves as the world coordinate. There is also a mobile Cartesian reference frame {P}, for which the origin is at the center P x p0 , y p0 , z p0 of the moving platform. A i (i = 1, 2, 3) are joints that are located at the center of the parallel chain of each arm and there is a pair of universal joints at its two ends and B i (i = 1, 2, 3) denote the location of joints on the moving platform. Three arms connect the joint A i and the joint B i and firmly connect the fixed base to the moving platform. In Figure 3, q i (i = 1, 2, 3) denote the actuated joint variables that are defined as the length between the point O(x 0 , y 0 , z 0 ) and the points A i . The slider moves when the PH-LDR operates, so q i are not constants. q max is the maximum value for q i and are equal to the maximum stroke of the rod-less PAs. Table 1 shows the geometry parameters for the PH-LDR. In this study, q i are defined as the joint variables; and the Cartesian variable vector O P P = [x y z] T represents the coordinate of the point P x p , y p , z p with respect to the frame {O}. The coordinates of the vector O −−→ OA i with respect to the frame {O} are: For simplicity, the unit vector of P − − → PB i is also defined as the vector → u i because the orientation of {P} is identical to the orientation of {O} for this design. The coordinate of the vector P − − → PB i with respect to the frame {P} is: where r is the distance between P x p0 , y p0 , z p0 and B i (i = 1, 2, 3). In general, the superscript can be ignored if the vector is expressed with respect to the frame {O}. In the design of the PH-LDR, the origin of the frame {O} shares the (x − y) plane with the origin of the frame {P}, and the orientation of {P} has an orientaion that is identical to the orientation of {O},

Inverse Position Kinematics
For the inverse position kinematics, the actuated joint variables ( 1, 2,3) i q i  are solved for a specific position on the moving platform. In terms of Figure 3, the vector-loop closure equation for each arm is: Substituting (1) and (2) into (3) yields: In Equation (5), there are two solutions for each actuator. However, only the negative square root solution satisfies the current assembly for the mechanism in which the three actuators are inclined inward from top to bottom. Thus, the inverse position kinematic equation for ith arm is:

Parameters
Description Value r The length between P x p0 , y p0 , z p0 and B i (i = 1, 2, 3) 60 mm l The length of kinematic chain 670 mm q max The maximum of the actuated joint variables q i (i = 1, 2, 3)(the maximum stroke) 500 mm

Inverse Position Kinematics
For the inverse position kinematics, the actuated joint variables q i (i = 1, 2, 3) are solved for a specific position on the moving platform. In terms of Figure 3, the vector-loop closure equation for each arm is: Substituting (1) and (2) into (3) yields: where l is the length of the kinematic chain. From Equation (4), the solution for q i is: Appl. Sci. 2020, 10, 3526 6 of 21 In Equation (5), there are two solutions for each actuator. However, only the negative square root solution satisfies the current assembly for the mechanism in which the three actuators are inclined inward from top to bottom. Thus, the inverse position kinematic equation for ith arm is: Substituting (1), (2) and (7) into (4) yields: As O P P = x y z T is given, the set of solutions for q i is:

Forward Position Kinematics
For a set of the input actuated joint variables, q i (i = 1, 2, 3), the position of the moving platform, O P P = [x y z] T , is solved using the forward position kinemactics. In Equation (8), combining the second equation with the first equation yields: Additionally, combining the third equation with the first equation yields: For q i (i = 1, 2, 3), the solution from Equations (10) and (11) is: Additionally, Substituting Equations (12) and (13) into the first equation in Equation (9) yields: Equation (14) gives two solutions for z. z must be negative, because the moving platform must move under the fixed base, so the solution of z is: Appl. Sci. 2020, 10, 3526 The unique solution of the forward position kinematic for the PH-LDR is defined in Equations (12), (13), and (15).

Velocity Analysis
T is the velocity vector for the sliders on the three axes and vector v = [v x v y v z ] T is the velocity vector for the moving platform, the Jacobian matrix J p between v q and v is defined as [21,22]: Differentiating Equation (8) with time yields: Its matrix form is: Additionally, From Equations (16) and (18), we have: where The determinant of J p is calculated as: Obviously, singularities occur for the following conditions: Contrastingly, Appl. Sci. 2020, 10, 3526 8 of 21 Equations (24) and (25) identify the singularities at which z = 0 or q 1 = q 2 = q 3 = r. There is no singularity at z = 0, because it implies that the moving platform is on the fixed base. Additionally, the initial length r of the moving platform is set to 60mm so the singularity q 1 = q 2 = q 3 = r possibly occurs when the three input actuated joint variables q 1 , q 2 and q 3 reach a distance of 60mm at the same time. In our design, the length of the parallel chain of each arm for the PH-LDR is designed as 130mm and the joints A i are located at the center of the parallel chain of each arm, so the minimum value for the joint varibles q i is 65 mm. Consequently, the singularity q 1 = q 2 = q 3 = r certainly can be prevented for the operation of the PH-LDR.

Controller Design
The LDR is a complex MIMO system with high nonlinearities and coupling states. In general, one can design a trajectory tracking controller for the end-effector through a kinetic model or a dynamic model. Thus, nonliearities and time-varying characteristics increase when the PH-LDR carries heavyweight and high mass moment of interia, and it is so difficult to model the exact dynamics for overall control systems of the PH-LDR. Moreover, the complicated and nonlinear dynamic model make the PH-LDR auduous to design a controller to manipulate the trajectory tracking in real time, due to massive numerical computation to solve partial differential equations. The purpose of this study is to implement the PH-LDR based on the kinetic model with the intelligence based controllers, which has the self-adaptation ability to deal with system nonlinearities and parameter uncertainties.
The dynamic for the PH-LDR has two parts: a rod-less PAs dynamic and a LDR motion dynamic. The rod-less PAs dynamic model can be separated into a pneumatic servo valve dynamic model and a pneumatic cylinder dynamic model, both of which have a significant effect on the stability of the PH-LDR, due to the highly nonlinear behavior caused by the airflow and air compression. In addition, without loss of generality, the motion dynamics of the LDR can be modeled in terms of the motion equations to configure the LDR in relation to the joint torque that is exerted by the three rod-less PAs. The inverse dynamics control (IDC), which is the so-called the computed torque control [23], is a well-known control design technique for conventional LDRs. It refers to the LDR motion model and directly predicts the desired output forces or torques for actuators using a feedforward method. However, the derived PH-LDR dynamic in this study is too complicated to obtain an accurate mathematical representation, so that the IDC method is difficult to be implemented in the developed PH-LDR.
Hence, a hybrid controller is proposed, based on a simple and straightforward technique that accounts for nonlinearities and uncertainties of PH-LDR under environmental temperature changes. The designed hybrid controller consists of a fuzzy-PID controller and a NNPC. Figure 4 shows the overall control system for the PH-LDR, in which each of the three independent Fuzzy-PID + NNPC's regulates the position of a rod-less PA. The fuzzy-PID receives the feed-back measured position and sends a control voltage to the proportional servo valve to regulate the amount of air that flows into the rod-less PA. The NNPC compensates for turnaround delay due to variations in the target, ∆Target and ∆ 2 Target, and variations in the temperature of the environment.
The inverse position kinematic model calculates the reference positions (ρ 1 (t), ρ 2 (t), ρ 3 (t)) for the three axes using a specific 3-DOF trajectory coordinate (x r (t), y r (t), z r (t)) . Accordingly, the feedback errors and the variations in the target are determined as: Additionally, where T s is the sampling time in the control loop. The outputs of (x(t), y(t), z(t)) define the 3-DOF position of the moving platform. Appl. Sci. 2020, 10, x FOR PEER REVIEW 9 of 21

Fuzzy-PID Controller Design
In this paper, the fuzzy-PID controller uses a fuzzy inference system to determine the variations of the PID controller parameters, i.e., ( , , ) The knowledge base for a fuzzy controller consists of a rule base and membership functions. The triangle-type fuzzy membership functions for ( )   Table 2, 3, and 4, these were designed by experts. In this paper, 49 fuzzy rules are created to implement the proposed fuzzy-PID controller for

Fuzzy-PID Controller Design
In this paper, the fuzzy-PID controller uses a fuzzy inference system to determine the variations of the PID controller parameters, i.e., ∆k i P , ∆k i I , ∆k i D , for the PH-LDR. The feedback error and error change, e i , . e i , were adopted for yielding the online self-adjustment ∆k i P , ∆k i I and ∆k i D . The alterations in the control parameters k i P (t), k i I (t) and k i D (t) in each time step are shown as follows: The knowledge base for a fuzzy controller consists of a rule base and membership functions. The triangle-type fuzzy membership functions for e i (t) and . e i (t) are shown in Figure 5a, and the membership functions for k i P (t) are shown in Figure 5b. The fuzzy membership functions for k i I (t) and k i D (t) are shown in Figure 5c. The fuzzy subset is {NB, NM, NS, ZO, PS, PM, PB}, in which the NB, NM, NS, ZO, PS, PM, and PB stand for negative big, negative medium, negative small, zero, positive small, positive medium and positive big, respectively. These levels of each parameter were determined based on the experimental properties and characteristics of the PAs. In this paper, expert's experience and knowledge are used to build a rule base. There are 49 rules in the fuzzy logic controller for ∆k i P , ∆k i I , and ∆k i D , respectively, shown in Tables 2-4, these were designed by experts. In this paper, 49 fuzzy rules are created to implement the proposed fuzzy-PID controller for the PH-LDR. The outputs of the fuzzy PID can be obtained using the Mamdani fuzzy inference and the center of area (COA) defuzzification.
the PH-LDR. The outputs of the fuzzy PID can be obtained using the Mamdani fuzzy inference and the center of area (COA) defuzzification.

NNPC Design
To prevent leaks from the high-pressure gas transmission, the gap between the rod-less PA and the slider is sealed using rubber bands. However, during the operation of PH-LDR, static friction and hysteresis between the slider and the pneumatic cylinder frequently produce a turnaround delay. A pre-compensator that adjusts the control force to account for variations in the reference position ρ i (t) is commonly used to address turnaround delay. This paper thus formulates the pre-compensation control as a linear combination of ∆Target and its rate of change, ∆ 2 Target, which is expressed as: where K 1 and K 2 are constant gains, V o is an offset pneumatic control voltage for the proportional servo valve, and the rate of change, ∆ 2 Target, is calculated by: When the slider moves on the rod-less pneumatic cylinder, the degree of frictional resistance to forward and backward motion is different, so the pre-compensation control is designed as: where K i,F and K i,N denote the respective compensation gains for the slider in the forward and backward directions. When the slider moves forward, ∆Target(t) > 0 so V PC (t) must be greater than the median voltage. On the contrary, when the slider moves backward, ∆Target(t) < 0 so V PC (t) must be smaller than the median voltage. The range for the offset pneumatic control voltage in the pre-compensation control V PC (t) is V o,F > 5V and V o,N < 5V, since the proportional servo valve works within the interval [0V, 10V]. When the parameters are fine-tuned, the pre-compensation control in Equation (31) is effective at dealing with the turnaround delay of each rod-less PA. However, the static friction and hysteresis in the PH-LDR exhibit highly nonlinear and coupling behavior and their values vary with time and operating temperature, so it is difficult to design accurate pre-compensation control using only the method in Equation (31). This study uses an artificial neural network (ANN) for NNPC to estimate the pre-compensation control V PC (t), in which the inputs of the ANN are the environmental temperature, ∆Target(t) and ∆ 2 Target(t), and the output of the ANN is V NNPC (t). As shown in Figure 6, the structure of the ANN that is used for this study is a feedforward propagation network with three layers, in which the input layer has three neurons, the hidden layer has ten neurons and the output layer has one neuron. Each neuron in the ANN calculates the weighted sum of the inputs from the previous layer and then produces an output with a bias value, as written by the following: where neuron i,k is the ith input node in the kth layer, W i,j,k is the weight between the jth node in the (k − 1)th layer and the ith node in the kth layer, O i,k is the output of the ith node in the kth layer, and Wb i,k is the bias of the ith node in the kth layer. where , i k neuron is the ith input node in the kth layer, is the weight between the jth node in the (k − 1)th layer and the ith node in the kth layer,  The NNPC for each rod-less PA uses the Matlab neural network toolbox to train the ANN and then the weights and bias for encoding the input-output relationship are obtained from training samples. The root-mean-square error (RMSE) is chosen as the performance function. To accelerate training speed and to avoid the problem of local minimum that is common to conventional back propagation (BP) algorithms, Levenberg-Marquardt BP is used and the training procedure will automatically terminate when the performance of generation cannot be improved. In this study, for each rod-less PA, a NNPC with a specific topology is trained and the ANN that performs best through the test set is selected. During the trajectory tracking operation of PH-LDR, the pre-compensation control for each rod-less PA is then calculated from the associated NNPC, in order to counterbalance the uncertain static friction and hysteresis. In such a case, the position regulation of rod-less Pas, even under the temperature change, can be well performed with the combination of Fuzzy-PID controller and NNPC in the feedback loop. Figure 7 shows the test-rig for the PH-LDR. Two experiments are used here to determine the feasibility of the PH-LDR. Figure 5 and Tables 2-4 respectively show the fuzzy membership functions and the fuzzy rules for the Fuzzy-PID controller. The temperature in the laboratory where the testrig was located was between 24 °C and 28 °C. The data sets for target  and 2 target  for the NNPC were collected between 24 °C to 28 °C, with an interval of one degree. 50% of the collected data was used to train the ANN and 50% was used for testing. During the training phase, the learning rate was set to 0.01 and the momentum coefficient was set to 0.91. The training was terminated if RMSE% < 10%, or after 10,000 iterations. In the present study, when the RMSE between the target value of ( ) PC V t and the estimated value ( ) NNPC V t is less than 0.1, the training procedure of ANN will be terminated. To evaluate the tracking performance for the PH-LDR, the instantaneous overall error of the moving platform is defined as:  The NNPC for each rod-less PA uses the Matlab neural network toolbox to train the ANN and then the weights and bias for encoding the input-output relationship are obtained from training samples. The root-mean-square error (RMSE) is chosen as the performance function. To accelerate training speed and to avoid the problem of local minimum that is common to conventional back propagation (BP) algorithms, Levenberg-Marquardt BP is used and the training procedure will automatically terminate when the performance of generation cannot be improved. In this study, for each rod-less PA, a NNPC with a specific topology is trained and the ANN that performs best through the test set is selected. During the trajectory tracking operation of PH-LDR, the pre-compensation control for each rod-less PA is then calculated from the associated NNPC, in order to counterbalance the uncertain static friction and hysteresis. In such a case, the position regulation of rod-less Pas, even under the temperature change, can be well performed with the combination of Fuzzy-PID controller and NNPC in the feedback loop. Figure 7 shows the test-rig for the PH-LDR. Two experiments are used here to determine the feasibility of the PH-LDR. Figure 5 and Tables 2-4 respectively show the fuzzy membership functions and the fuzzy rules for the Fuzzy-PID controller. The temperature in the laboratory where the test-rig was located was between 24 • C and 28 • C. The data sets for ∆target and ∆ 2 target for the NNPC were collected between 24 • C to 28 • C, with an interval of one degree. 50% of the collected data was used to train the ANN and 50% was used for testing. During the training phase, the learning rate was set to 0.01 and the momentum coefficient was set to 0.91. The training was terminated if RMSE% < 10%, or after 10,000 iterations. In the present study, when the RMSE between the target value of V PC (t) and the estimated value V NNPC (t) is less than 0.1, the training procedure of ANN will be terminated. To evaluate the tracking performance for the PH-LDR, the instantaneous overall error of the moving platform is defined as:

Experiments
where e x (t), e y (t), e z (t) is defined as (x d (t) − x(t)), (y d (t) − y(t)), (z d (t) − z(t)) . The position (x(t), y(t), z(t)) is calculated by the forward position kinematics, as q i (i = 1, 2, 3) are given. e 1 (t), e 2 (t) and e 3 (t) respectively represent the tracking errors for the 1st-axis, the 2nd-axis and the 3rd-axis. It is proven in Section 3 that there are no singularities on the workspace for the PH-LDR. In this paper, the Monte Carlo method was used for the workspace analysis of the PH-LDR. The Monte Carlo method is also called the extreme limit boundary searching method, which is conducted based on the inverse solution of the kinematic location of the PH-LDR. Given all the possible workspaces of the PH-LDR, a large number of random points will be generated, and every point will be judged whether to be in the workspace or not. If it meets the constraint condition, it belongs to the workspace; if not, it shall be rejected. All the feasible points will constitute the workspace of the PH-LDR. The one from the qualified to the unqualified is just the boundary point of the workspace, and the line constituted is the boundary line. The structural parameters of the PH-LDR are shown in Table 1. The search step length of each variable was set to be 2mm. Through the Equation (9) and employment of MATLAB, the 3D stereoscopic graph of the workspace of this PH-LDR was plotted, as shown in Figure 8. The detailed steps of the Monte Carlo method for the workspace can refer to [24].
Appl. Sci. 2020, 10, x FOR PEER REVIEW 13 of 21 axis. It is proven in Section 3 that there are no singularities on the workspace for the PH-LDR. In this paper, the Monte Carlo method was used for the workspace analysis of the PH-LDR. The Monte Carlo method is also called the extreme limit boundary searching method, which is conducted based on the inverse solution of the kinematic location of the PH-LDR. Given all the possible workspaces of the PH-LDR, a large number of random points will be generated, and every point will be judged whether to be in the workspace or not. If it meets the constraint condition, it belongs to the workspace; if not, it shall be rejected. All the feasible points will constitute the workspace of the PH-LDR. The one from the qualified to the unqualified is just the boundary point of the workspace, and the line constituted is the boundary line. The structural parameters of the PH-LDR are shown in Table 1. The search step length of each variable was set to be 2mm. Through the Equation (9) and employment of MATLAB, the 3D stereoscopic graph of the workspace of this PH-LDR was plotted, as shown in Figure 8. The detailed steps of the Monte Carlo method for the workspace can refer to [24].   axis. It is proven in Section 3 that there are no singularities on the workspace for the PH-LDR. In this paper, the Monte Carlo method was used for the workspace analysis of the PH-LDR. The Monte Carlo method is also called the extreme limit boundary searching method, which is conducted based on the inverse solution of the kinematic location of the PH-LDR. Given all the possible workspaces of the PH-LDR, a large number of random points will be generated, and every point will be judged whether to be in the workspace or not. If it meets the constraint condition, it belongs to the workspace; if not, it shall be rejected. All the feasible points will constitute the workspace of the PH-LDR. The one from the qualified to the unqualified is just the boundary point of the workspace, and the line constituted is the boundary line. The structural parameters of the PH-LDR are shown in Table 1. The search step length of each variable was set to be 2mm. Through the Equation (9) and employment of MATLAB, the 3D stereoscopic graph of the workspace of this PH-LDR was plotted, as shown in Figure 8. The detailed steps of the Monte Carlo method for the workspace can refer to [24].

Experiment 1-Star Loop Motion Trajectory (a Two-Dimensional Motion)
In Experiment 1, the PH-LDR tracks the plane trajectory of a star loop motion trajectory. The tracking objective is a star loop, which is composed of five oblique straight lines. The profile and the direction of movement for the trajectory are shown in Figure 9. The formulation is: (35) Appl. Sci. 2020, 10, x FOR PEER REVIEW 14 of 21

Experiment 1-Star Loop Motion Trajectory (a Two-Dimensional Motion)
In Experiment 1, the PH-LDR tracks the plane trajectory of a star loop motion trajectory. The tracking objective is a star loop, which is composed of five oblique straight lines. The profile and the direction of movement for the trajectory are shown in Figure 9. The formulation is: time that is required to travel along each of the segments 1 to 5 is 2 s. Figure 10 shows the desired star loop motion trajectory and the movement of the moving platform and the overall error 3 ( ) D e t . Figure 11 shows the desired trajectories and the system responses for the three rod-less PAs. Figure 12 shows the control voltages of the proportional servo valve for the 1st-axis, the 2nd-axis, and the 3rd-axis rod-less PA when the end-effector moves along the star loop motion trajectory. The proportional servo valve can adjust a spool according to inputvoltage to change compressed air that flows into the chambers of a rod-less PA. Initially, the inputvoltage is set to 5 V to remain the piston at the center of the rod-less PA in the experiments. Figure 13 shows the tracking errors of the 1st-axis, the 2nd-axis and the 3rd-axis rod-less PAs. The tracking error for the 1st-axis is within ±0.6 cm; the tracking error for the 2nd-axis is within ±0.58 cm; the tracking error for the 3rd-axis is within ±0.6 cm; and the maximum of 3 ( ) D e t for the star loop motion is less than 0.65 cm. There is a large tracking error in the start loop trajectory when the command position is suddenly changed; otherwise, the tracking error is very small. The results show good contour tracking for the star loop motion trajectory. In Equation (35), S is the position increment, x i is the x-axis trajectory, and y i is the y-axis trajectory of the moving platform. The actuated joint trajectories for the PH-LDR are computed by the inverse kinematics. To prevent the sliders and the shock absorber from touching at the ends of the 1st-axis, the 2nd-axis and the 3rd-axis rod-less PAs stroke during the trajectory tracking, the center position of the star motion trajectory is shifted by a bias that is larger than the amplitude of the trajectory. Initially, the moving platform moves along the fifth-order path in 3 s and the moving platform moves from the zero position to the starting point P 1 = (10, 0, −60cm) for the star motion trajectory. For segment 1 to segment 5, the platform moves along a star loop path with edges of length S = 19cm through five positioning points at the vertices of the star loop. The trajectory comprises six segments, which define five positioning points at the vertices of the stat: P 1 = (10, 0, −60cm), P 2 = (−8.1, 5.9, −60cm), P 3 = (3.1 − 9.5, −60cm), P 4 = (3.1, 9.5, −60cm) and P 5 = (−8.1, −5.9, −60cm). The time that is required to travel along each of the segments 1 to 5 is 2 s. Figure 10 shows the desired star loop motion trajectory and the movement of the moving platform and the overall error e 3D (t). Figure 11 shows the desired trajectories and the system responses for the three rod-less PAs. Figure 12 shows the control voltages of the proportional servo valve for the 1st-axis, the 2nd-axis, and the 3rd-axis rod-less PA when the end-effector moves along the star loop motion trajectory. The proportional servo valve can adjust a spool according to input-voltage to change compressed air that flows into the chambers of a rod-less PA. Initially, the input-voltage is set to 5 V to remain the piston at the center of the rod-less PA in the experiments. Figure 13 shows the tracking errors of the 1st-axis, the 2nd-axis and the 3rd-axis rod-less PAs. The tracking error for the 1st-axis is within ±0.6 cm; the tracking error for the 2nd-axis is within ±0.58 cm; the tracking error for the 3rd-axis is within ±0.6 cm; and the maximum of e 3D (t) for the star loop motion is less than 0.65 cm. There is a large tracking error in the start loop trajectory when the command position is suddenly changed; otherwise, the tracking error is very small. The results show good contour tracking for the star loop motion trajectory.

Experiment 2-Helical Circular Motion Trajectory
In Experiment 2, a helical circular motion trajectory is used to conduct space trajectory tracking control for the PH-LDR. The tracking objective is to force the moving platform along the helical circular motion trajectory. To prevent the slider and the shock from touching at the ends of the stroke of the 1st-axis, the 2nd-axis and the 3rd-axis rod-less PAs during trajectory tracking, the center position of the specified helical circular trajectory is shifted by a bias that is larger than the amplitude of the trajectory. The overall trajectory has two segments: a fifth-order path that takes the slider from the zero position to the starting point of the helical circular motion trajectory, and a segment wherein the platform moves along a helical circular motion trajectory. The helical circular motion trajectory is shown in Figure 14.
The trajectory is formulated as [25]: where θ i (t) = θ i (t − 1) + ∆θ, and ∆θ, ∆z, γ, x(t), y(t), z(t) are the increment of the angle, the increment of the z-axis, the radius, the x-axis position, the y-axis position, and the z-axis position of the moving platform, respectively. The helical circular tracking motion has a trajectory that is described by Equation (36), where the radius γ = 10 cm, θ i (t) are within the interval of [0 10π], ∆θ = π/300, ∆z = 1/50 mm and z(t) are from −65 cm to −60 cm. The desired path and the system response for the moving platform for the helical circular trajectory are shown in Figure 15. The maximum of e 3D (t) for the overall helical circular motion is less than 0.61 cm. Figure 16 shows the desired trajectories and the system responses for the three rod-less PAs. Figure 17 shows the control voltages of the proportional servo valve for the 1st-axis, the 2nd-axis, and the 3rd-axis rod-less PA when the end-effector moves along the star loop motion trajectory. The proportional servo valve can adjust a spool according to input-voltage to change compressed air that flows into the chambers of a rod-less PA. Initially, the input-voltage is set to 5 V to remain the piston at the center of the rod-less PA in the experiments. Figure 18 shows the tracking errors of the 1st-axis, 2nd-axis and 3rd-axis for the rod-less PAs when the platform moves along the helical circular trajectory, where the tracking error for the 1st-axis is within ±0.6 cm, the tracking error for the 2nd-axis is within ±0.7 cm, and the tracking error for the 3rd-axis is within ±0.7 cm. The tracking error in the 3-DOF motion contour tracking for the overall helical circular motion is small. Figures 15-18 show that the proposed rod-less PAs based PH-LDR is effective and accurate.
1 / 50 z mm   and ( ) z t are from −65 cm to −60 cm. The desired path and the system response for the moving platform for the helical circular trajectory are shown in Figure 15. The maximum of 3 ( ) D e t for the overall helical circular motion is less than 0.61 cm. Figure 16 shows the desired trajectories and the system responses for the three rod-less PAs. Figure 17 shows the control voltages of the proportional servo valve for the 1st-axis, the 2nd-axis, and the 3rd-axis rod-less PA when the end-effector moves along the star loop motion trajectory. The proportional servo valve can adjust a spool according to input-voltage to change compressed air that flows into the chambers of a rod-less PA. Initially, the input-voltage is set to 5 V to remain the piston at the center of the rod-less PA in the experiments. Figure 18 shows the tracking errors of the 1st-axis, 2nd-axis and 3rd-axis for the rod-less PAs when the platform moves along the helical circular trajectory, where the tracking error for the 1st-axis is within ±0.6 cm, the tracking error for the 2nd-axis is within ±0.7 cm, and the tracking error for the 3rd-axis is within ±0.7 cm. The tracking error in the 3-DOF motion contour tracking for the overall helical circular motion is small. Figures 15-18 show that the proposed rodless PAs based PH-LDR is effective and accurate.

Conclusions
This study develops a PH-LDR that consists of three rod-less PAs in a horizontal structure, so that a larger workspace can be yielded for the end-effector of LDR to achieve 3-DOF trajectory tracking. The parallel mechanism of the PH-LDR is analyzed based on a geometric method to determine the closed-form solutions of the kinematic chains. Moreover, a full-scale test-rig is experimentally set up in order to verify the feasibility of the proposed PH-LDR. For realizing the trajectory tracking control, a hybrid controller consisting of fuzzy-PID controller and NNPC is proposed, to overcome the nonlinearities and uncertainties resulted from the air compressibility and temperature change in the motion control of the pneumatic actuating system. From the experimental results, the design and construction of the PH-LDR system is verified through complex 3-DOF motion tracking control, including a star loop trajectory and a helical circular trajectory, on the developed test rig. In addition, the proposed architecture of LDR with the control methodologies can also be easily extended to the other 3D parallel mechanism robot, with modification of the workspace and the drive system.

Patents
This section is not mandatory, but may be added if there are patents resulting from the work reported in this manuscript.