Abstract

The determination of kinematic parameters for a parallel robotic system (PRS) is an important and a critical phase to maximize reachable workspace while avoiding singular configurations. Stewart Platform (SP) mechanism is one of the widely known PRS and it is used to demonstrate the proposed technique. In the related literature, GCI (Global Condition Index) and LCI (Local Condition Index) are the commonly used performance indexes which give a measure about the dexterity of a mechanism. In this work, Sequential Quadratic Programming (SQP) method is used to optimize kinematic parameters of a 6dof 3x3 UPU SP in order to reach maximum workspace satisfying small condition numbers. The radius of mobile and base platforms and the lengths of the legs used in the platform are chosen as kinematic parameters to be optimized in a multiobjective optimization problem. Optimization is performed at different stages and the number of optimized kinematic parameters is increased at each level. In conclusion, optimizing selected kinematic parameters at once by using SQP technique presents the best results for the PRS.

1. Introduction

In 1947, Gough [1] established the basic principles of a mechanism with a closed-loop kinematic structure that allows positioning and orientation of a mobile platform to test tire wear and tear. In 1965, Stewart offered a very similar structure to Gough’s mechanism to be used as an aircraft simulator. Ironically, the Gough platform is most often known as Stewart platform [1]. In recent years, there have been many researches about PRS whose principles are originated from Gough. They have numerous advantages compared to the serial robotic systems such as high load carrying capacity and precise positioning ability. A PRS has a mobile platform connected to a fixed base by parallel legs instead of serial links.

Besides their advantages compared to serial manipulators, PRS have a small workspace and require a more complex mechanical design. Controlling the motion of a parallel mechanism is a challenging task because of its complex kinematic structure and variable dynamics through its workspace. To overcome this issue, many researchers tried to find alternative parallel structures like redundant PRS [24]. In recent years, many different parallel robotic structures have been investigated for various applications.

The performance of a PRS is mainly affected by its geometrical parameters. Determining geometrical parameters is called geometric synthesis. One of the most important measures for kinematic performance is called dexterity which describes the ability of a mechanism to move and apply forces in arbitrary directions as easily as possible. However, the most suitable index to define the accuracy and rigidity performance is the condition number, which LCI (Local Condition Index) and GCI (Global Condition Index) are related to. LCI is defined as the ratio of the maximum and minimum singular values of the Jacobian matrix of a PRS [5]. It varies between 1 and infinity [6, 7]. GCI index is a measure of dexterity and the average of LCI calculated for entire workspace.

Dexterity and workspace are the most important criteria for optimal designs of parallel manipulators. There are many researches investigating the effects of geometric parameters on dexterity, workspace, and singularity. Some of them are about maximizing well-conditioned workspace [8] using different optimization methods. Monte Carlo method [8, 9], multiobjective/genetic algorithm [1012], controlled random search (CRS) [13, 14], analytic algorithm [15], and neurogenetic algorithm are mostly used ones [16]. G. Abbasnejad offered a practical algorithm to find the singularity free cylindrical workspace of parallel manipulators using particle swarm optimization method (PSO) [17, 18]. In the literature, there are other performance indexes such as GSI (Global Stiffness Index) and VTF (Velocity Transformation Factor) [19]. Therefore, some researchers optimized PRS for multiperformance indexes which are called MOO (Multiobjective Optimization) [2022]. These optimization methods are used for investigating the best geometric parameters for structures of parallel manipulators. On the other hand some researchers compared convergence performances of optimization methods [23].

In this research, we used SQP method which is Lagrangian based [24] for the optimization of kinematic parameters to reach the maximum workspace for a PRS. This paper is organized as follows. In the following section, the mechanism of a PRS and its geometric parameters and Jacobian matrix formulation for 3x3 SP are briefly described. In Section 3, our performance indexes and workspace assumption are described. Section 4 illustrates the SQP and optimization algorithm. In the last section, we discussed the results of the SQP optimization procedure.

2. Mechanism of the PRS

A PRS is constructed from parallel connected legs between mobile and base platforms. There are two major classification criteria of a PRS which are based on (a) joint types of connections of legs to the mobile and base platform and (b) the number of legs. The array of the types of joints are universal, prismatic, and universal, respectively. So it is called UPU. Specifically, linearly driven six legs which are connected to a mobile and a base platform cause a full spatial motion of 3-axis translational and 3-axis rotational motion. This PRS is called Stewart Platform [20, 25, 26] and classified by the number of leg connections on upper and lower platforms. A 3x3 SP includes 3 connections on lower and upper side of the legs. Each two of six legs share cocentered joints on each side of the legs. “Figure 1” illustrates the general structure of SP that 3x3 special configuration is used for the current research.

An SP has 48 design parameters affecting its workspace and kinematic behavior [27]. These parameters are spatial positions of upper and lower joints (2x18 parameters) and upper/lower leg lengths (2x6 parameters). By placing connection points on a circle, these would reduce to upper/lower circle diameters/radiuses passing through upper/lower joint centers, rp and rb. In this research, optimization parameters are these two radiuses and length of the legs, .

2.1. Geometric Parameters of PRS

The lengths of the legs (), which change according to the desired position and orientation of the mobile platform, can be calculated using the inverse kinematics equations, x (x, y, z, α, β, γ). In this case, there are six unknowns, . “Figure 1” helps us understand the geometry of SP. Mobile and base coordinate frames, and , are defined. They each are linked to the mobile platform and the base. The origin of frame is located at the center of a circle passing through the center of joints between the mobile platform and the legs. The axis is pointing upward, and the axis is the bisecting line of the angle (). The origin of frame is located at the center of a circle passing through the center of joints between the base and the legs. The axis is pointing upward, and the axis is the bisecting line of the angle .

These two coordinate frames lead us to determine the position of the links attached to the base at (B1…, B6) and the mobile platform at (P1.., P6).

The angle (, ) is denoted as . Angles (, ), (, ), and (, ) are each equal to 120°. Similarly, angle (, ) is denoted by . Angles (, ), (, ), and (, ) are, also, equal to 120°.

Now, the length of the legs () can be calculated forThe angles and are denoted as and , respectively. Therefore, we have “Eq.” (5) represents translational position of with respect to ,“Eq.” (6) is representing the rotation of mobile coordinate frame with respect to base frame .“Eq.” (7) describes joint positions of the mobile platform with respect to .“Eq.” (8) describes joint positions of the base platform with respect to .

, , , , , where is the rotation of with respect to about . is the rotation of with respect to about . is the rotation of with respect to about .

Consequently,The length of the legs can be calculated as follows:Using (4), we get length of the legs as

2.2. Jacobian Matrix

Jacobian matrix of a robotic system represents the linear transformation between the speeds in joint space, , and world space, , of a robotic manipulator. Jacobian matrix provides the velocity ratio from the joint space to the space of the end-effector [26].

For a parallel manipulator, the length of the legs () as a vector is given asThe position of the end-effector isSo, the relation equation between the position of the end-effector (x) and the length of the legs () can be represented asBy differentiating (15), we obtainConsequently is the Jacobian matrix of the mechanism, which is used in calculating the performance indexes of the PRS, and is utilized in this research as

3. Performance Indexes and Workspace

3.1. Dexterity, GCI, and LCI

Dexterity is an indicator of the ability of a parallel manipulator about how rates of joints transform to the rates of the mobile platform and vice versa. It depends on the position and the structure of the parallel manipulator and can be determined by using the condition number of its Jacobian matrix aswhere Cond (J) is the condition number of Jacobian matrix and is defined as represents the norm of the related vector or matrix.

While calculating the condition number of a Jacobian matrix, different types of norms are used in the literature, and the most used ones are Frobenius and Spectral norms [1]. The difference between these norms is about their calculation methods. If the Frobenius norm is considered, then we can writeIn this case, the Frobenius norm is defined as the extracting roots of the quadratic sum for each element of the Jacobian matrix. GCI (Global Condition Index) is defined all over the workspace as the average of Cond(J) = LCI (Local Condition Index).

If the spectral norm is introduced, the dexterity indexes can be described aswhere and represent the maximum and minimum singular values of the Jacobian matrix (J), respectively. The values of Cond (J), which are directly related to the singular values of the Jacobian matrix, are between one and positive infinity. All singular values will be the same in the Jacobian matrix, and the manipulator is isotropic when Cond (J) is equal to 1.

Otherwise, when Cond (J) is potentially equal to positive infinity, then the Jacobian matrix is singular. So, when the condition number is considered, the value of the condition number should be minimum for a better dexterity.

The condition number is used as an index for measuring the accuracy/dexterity of a robotic system, and also it possesses information of whether a position is close to singularity or not. However, it is not entirely possible to define a mathematical proximity to singularity for a PRS whose degree of freedom is a mix of translational and rotational motion. This is because the use of condition number is just an index and yet it has an advantage as a single number which is used to describe the overall kinematic behavior of a PRS.

Because of the complex definition of the condition number, we cannot calculate its analytical form as a function of the position parameters except for very simple robots. However, robust linear algebra tools can calculate it numerically for a given spatial position. Parallel manipulators that have both translational and rotational dof have a major drawback about condition number. Because the unit of translational motion is different from the rotational motion, the Jacobian matrix is not homogenous [28]. Thus, there exist some researches about the normalization of Jacobian matrices which are partly rotational and translational. Some of them defined their own methods for normalization and some of them divided Jacobian matrix by the degree of the matrix [2932].

3.2. Uniformity

Uniformity is another global performance index for PRS. It is an indicator of the uniform dexterous workspace. Calculation of uniformity index is shown in equation (25). Uniformity index value changes between 1 and infinity. The index value of the best uniform dexterous workspace is equal to 1.

3.3. Workspace

The number of possible mechanical configurations of serial robotic systems is less than that of parallel robotic system configurations. Thus, comparison of workspace of serial robotic systems is applicable. For example, one can easily assess that RRR robotic system has greater workspace capacity than a cartesian configuration because of its angular workspace capacity. However, because of the vast number of mechanical configuration possibilities of parallel robotic systems, no such comparison could be performed easily [24].

One can manage to compare 6dof parallel robotic systems that are abbreviated as SP in this study. These classes of full-motion-capable spatial mechanisms are classified as 3x3, 6x3, and 6x6 based on their number of upper and lower platform connection points. So, 3x3 has triangular base and mobile platform, while 6x3 configuration has hexagonal base and triangular mobile platform configuration. Lastly, 6x6 has hexagonal base and mobile platform. Considering the first, the second, and the third configuration of similar dimensions and with the same minimal and maximal leg lengths, the third has the largest workspace, followed by the second and then the first. It must be noted that the workspace volume of the 6x6 structure is always approximately 30% greater than that of the 6x3 mechanical structure. Similarly the 6x6 configuration volume is 70% greater than that of the 3x3 workspace. Thus the 6x3 volume is approximately 25% greater than that of the 3x3 configuration workspace [1].

Additionally, the effect of geometrical parameters on parallel robotic systems is much greater than that on serial robotic systems. For instance, external stiffness over a given workspace of a 6x6 Gough Platform can differ by 700% for a change of only 10% of the platform radius, [24].

Therefore, in order to extend the reachable workspace of a 3x3 SP mechanism, objective workspace for optimization and performance evaluation is described as follows:(i) = Changes from -0.06 to +0.06 m in 0.04 m regular intervals(ii) = Changes from -0.06 to +0.06 m in 0.04 m regular intervals(iii) = Changes from +0.3 to +0.4 m in 0.02 m regular intervals(iv), , = Changes from -5 to +5 in 2.5° regular intervals

Therefore, total 12000 spatial, 6dof, coordinates need to be satisfied after the optimization performance. First, according to defined workspace, we calculated LCI in Frobenius norm for 12000 points in the workspace and then calculated GCI. “Figure 2” represents the volume of the workspace for 3x3 SP and the inner workspace inside of it to be optimized. We used different resolution, “Δa”, values for each axis of the workspace as mentioned above. For example, Δa is 0.04 m for linear motion along y-axis.

4. Optimization by SQP

The primary purpose of the optimization is to increase dexterity and usable workspace of the platform by changing the design parameters. Besides various optimization techniques, we used SQP method for the optimization of design parameters , , and . Min and max values of the objective design parameters are set before the optimization procedure. Another constraint is Z (LCI of workspace point). We calculated all Z (LCI) for each point in the predefined workspace shown in “Figure 2.”

4.1. Cost Function Assignments

In this research, Matlab Optimization Toolbox (fmincon) is used to perform SQP optimization technique for the PRS. The optimization procedure is accomplished in three stages for maximum dexterity with reachable workspace:(1)One-variable optimization: is intended to be optimized(2)Two-variable optimization: and are intended to be optimized(3)Three-variable optimization: , , are optimized all together

For all optimization approaches, a sub cost function is defined as follows:(i)Sub cost function for one-variable optimization is(ii)Sub cost function for two-variable optimization is(iii)Sub cost function for three-variable optimization is(iv)Cost function for all cases is

We define cost functions as described above because minimum g(x) values are investigated, which occur from the sum of minimum sub cost function, f(x), values. It means that minimum f(x) values are to be found indeed. While searching minimum f(x) values, two criteria are present. One of them is the need for minimum Z and minimum and the second criterion is the reachable workspace. Since Z never becomes smaller than 1, the algorithm is constructed accordingly. If the targeted position in the workspace is not reachable, large value of Z is assigned not to calculate real Z value in order to eliminate the unreachable point for evaluation. Thus, unnecessary calculations finding minimum f(x) for that corresponding position are avoided. Searching for the minimum of two variables, and , which makes the sub cost function f(x) minimum, we multiply these variables defining the sub cost function . Thus, while minimizing the cost function, we optimize variables. We constrained Z (condition number) with a maximum value of 1000 as an indicator for dexterity. We limited leg lengths and checked if the coordinates in the workspace are reachable or not. Searching for the minimum of three variables , , and , which make the sub cost function f(x) minimum, we multiply , with an average of and (Z-1) defining the sub cost function . We optimized parameters for maximum reachable dexterous workspace.

4.2. One-, Two-, and Three-Variable Optimization Algorithm

Figure 3 represents the optimization procedure sequence for one variable in which we have to fix other variables and define their min and max values to constrain them. Secondly, we limit the leg lengths to check if the points in the workspace are reachable or not. Later on, we define the sub cost function according to one variable or for two variables and create a limit value for selected . After that, we check if the workspace points are reachable or not, and we eliminate the unreachable ones. Then, cost function is defined and condition number for each workspace point is calculated and an optimum value for is investigated which makes the sum of (or ) minimum. A similar algorithm with a little difference is performed for two-variable optimization, and . When optimizing for two variables, first step defined in “Figure 3” is skipped and instead of that step, we define min-max values for as a second stage, and we use sub cost function set for two variables instead of at the fifth step. For three-variable optimization, we constrained all variables defining their minimum and maximum values. As a sub cost function, we used .

5. Optimization Results

3x3 SP mechanism mobile platform diameter is analysed for constant base diameter utilizing SQP technique. Mobile platform/base platform diameter optimization for better workspace reach and dexterity is performed. Condition number variation for different mobile platform diameters in selected range, minimum of condition number for base/mobile platform diameter change, and the number of reachable points for different base/mobile platform diameter combinations are discussed.

In Figure 4, horizontal axis refers to scanned values, and vertical axis refers to minimum condition number which we get from calculations without homogenization of Jacobian matrix. We get low values of while is increasing. We are searching optimum values for most efficient workspace volume which has reachable points with low condition numbers. As a result, high mobile platform diameter refers to calculated minimum condition numbers.

In Figure 5, horizontal axis refers to scanned values. Left vertical axis refers to scanned values and right vertical axis refers to the mapping of condition numbers which we calculated in Frobenius norm for predefined reachable workspace. As seen in Figure 5, minimum condition numbers occur with high values of and which mean dexterity is getting better while they are increasing. However, it is obvious from the plot that comparable decrease exists of the mobile platform diameter with respect to base diameter to obtain minimal condition numbers. Also, the rate of this decrease is getting smaller when the base diameter is a lot higher than mobile platform diameter. The horizontal tendency of platform diameter is apparent for higher base diameters.

In Figure 6, horizontal axis refers to scanned values. Left vertical axis refers to different values, and the right vertical axis refers to the mapping of reachable points in the workspace which have condition numbers limited with an upper bound of 1000. As indicated in Figure 6, while values decrease, values increase, and numbers of reachable points elevate. The upper bound of Z value is equal to 1000. Thus, the numbers of reachable points are high for low mobile platform diameters and high base diameters.

During optimization process of mobile platform diameter, different initial values of with various constraints are used, as shown in Table 1. While = 0.175 [m], minimum limit of is kept 0.050 [m] and leg lengths are scanned between 0.300 [m] and 0.450 [m]. As a result, the optimum value of is calculated to be 0.075 [m]. When minimum limit for is 0.080 [m], optimum value is found as = 0.088 [m]. As mentioned above, different initial values of are used for optimization. Consequently, when the constraint value of changes, the optimum value of changes. This is because, in different workspace regions, we have different local minimums.

If the optimization procedure in Figure 8 is analyzed thoroughly, while the lowest constraint of is 0.050 and 0.070 [m], we selected optimization scanned range between 0.050 and 0.180 [m]. Then, we obtain an optimum value of 0.075 [m]. When the min constraint value of is elevated to 0.080 [m], the first local minimum of 0.075 [m] is passed. Optimized value of shifts from 0.075 [m] to 0.088 [m]. Out of nine different optimization combinations, one can understand that we have two local minimums which are 0.075 [m] and 0.088 [m]. Then, combinations that include the initial value of 0.090 [m] and larger will result in an optimized value of start point of the scanned range, 0.090 [m], 0.100 [m], 0.110 [m], etc.

Consequently, we found two local minimum values (0.075 [m] and 0.088 [m]) for having maximum reachable points which are limited to upper bound condition number, Z, of 1000. Optimization results for two variables and and their cost function values with different initial values of the optimization procedure are as shown in Table 2. For example, while the initial value of is equal to 0.125 [m] and an initial value of is equal to 0.070 [m], cost function value is equal to 4.8xE7, optimum value of is equal to 0.074 [m], and optimum value of is equal to 0.12 [m]. Optimum values of diameters and depend on constraints.

While generating Table 2, constraints are rb: 0.125 [m] < rb < 0.175 [m] and constraints for ; 0.070 [m]  [m]. Low mobile platform diameters with high base diameters seemed to result comparably in small cost function values, which means better dexterity in motion. In Figure 7, minimum cost function value is 4.8 x 107 as indicated in bold areas. In those areas, cost function values are minimum with optimized and values as shown and are compatible with Table 2. While optimized has minimum values, cost function values are minimum. As mentioned previously, we are searching for minimum cost function value for better dexterity. In Table 2, there is not one optimum solution for the minimum cost function value. For example, , = 0.074 [m], 0.129 [m] has minimum cost function value of (4.8xE7) and also , = 0.077 [m], 0.125 [m] has minimum cost function value of (4.8xE7).

For better accuracy, we increase the resolution in solution area in Table 3. For the purpose of interpreting the results, we used cost function values. However, the condition number has large spectrum throughout the workspace. Therefore, comparing the average condition numbers for different design parameters with different constraints is meaningless. As a result, calculating the cost function (sum of condition numbers multiplied by rb and rp) and considering the resulting value is more efficient while searching for the best optimum design parameters.

Optimum values are in Table 3. Minimum cost function value is 3.2xE7. Solution set with minimum cost function value occurs for , = 0.072 [m], 0.127 [m]; , = 0.071 [m], 0.127 [m]; , = 0.071 [m], 0.128 [m].

Optimization results for three variables , , , and their cost function values with different initial values of the optimization are as shown in Tables 4 and 6. For example, while the initial value of is equal to 0.125 [m], initial value of is equal to 0.070 [m], and initial value of is equal to 0.300 [m], cost function value is equal to 4.8xE7, optimum value of is equal to 0.076 [m], is equal to 0.125 [m], and is equal to 0.300 [m]. Optimum values of parameters also depend on constraints. While generating Table 4, our constraints are as in Table 5.

In Table 4, there is not only one optimum solution for the minimum cost function value. For example, , , = 0.076 [m], 0.125 [m], 0.300 [m] has minimum cost function value of (4.8xE7) and also , , = 0.077 [m], 0.125 [m], 0.388 [m] has minimum cost function value of (4.8xE7). For better accuracy, we zoom into the solution area in Table 5. So, we generated Table 6. As shown in Table 6, minimum cost function value is found as 3.2xE7. According to Table 6, optimum parameters are , , = 0.071 [m], 0.127 [m], 0.350 [m] and , , = 0.072 [m], 0.127 [m], 0.374 [m]. Optimization results change according to optimization parameters which are initial values and convergence tolerances. Hence, we controlled our optimization results with mathematical calculations and compared these values with calculated values.

For three-parameter optimization, we used as a sub cost function and obtained results in Tables 4 and 6 using the optimization algorithm as shown in Figure 8.

For Table 7, we used optimized values from Table 3, and we generated Table 7 using inverse kinematic equations.

Table 7 presents Lmax, Lmin, (uniformity index), and reachable points for optimum radius values of mobile and base platforms. According to these calculated values using optimized parameters from Table 7, GCI (Global Condition Index) and LCI (Local Condition Index) values are computed for each solution. Figure 9 represents LCI values for the first optimized variables. As shown in Figure 9, there are nearly 400 unreachable points in the workspace. LCI values of reachable points vary between nearly 8 and 11. The average of LCI values is defined as GCI.

When is equal to 0.127 [m], is equal to 0.071 [m], and is between 0.300 and 0.450 [m], LCI values are as shown in Figure 9 and GCI is equal to 9.62.

When is equal to 0.127 [m], is equal to 0.072 [m], and is between 0.300 and 0.450 [m], LCI values are as shown in Figure 10 and GCI is equal to 9.49.

When is equal to 0.128 [m], is equal to 0.071 [m], and is between 0.300 and 0.450 [m], LCI values are as shown in Figure 11 and GCI is equal to 9.57.

As shown in Figures 9, 10, and 11, LCI and GCI values are very similar, because two parameters and differ only 0.001 [m]. Distribution of LCI values in graphics is different, but it does not affect GCI noticeably. As shown in Table 7, we can change constraints of and reach all workspace. One can redefine leg lengths of the PRS. If we compare Tables 6 and 7, calculated solutions are different but very similar to each other because initial values affect optimization procedure and convergence tolerances of optimization. All LCI in Figures 9, 10, and 11 are calculated from the homogenous Jacobian matrix and these calculated values are divided by the dimension of Jacobian matrix for homogenization.

In detail, first solution is the best for the amount of the reachable points. Second one is the best solution according to GCI. And the third one is the best for the uniformity.

6. Conclusion

Determining the values of geometrical parameters is very important for efficiency of a PRS design. In literature, there are many researches for the optimization of kinematics of a PRS. In these works, various optimization methods have been used. In this study, we optimized 3x3 SP mechanism using SQP optimization method for predefined workspace and checked the efficiency of the optimization for three parameters which are diameter of mobile and base platform together with lengths of the legs. It is shown that the optimization of multi-parameters is not easy to formulate concerning cost functions and results must be evaluated by calculating the dexterity indexes and the amount of reachable points. Firstly, we optimized and then, we optimized and simultaneously. Lastly, , , and are optimized at the same time. During optimization, tolerances, procedure steps, and applied/defined cost functions affect the convergence efficiency. Results show that initial values for each parameter at the beginning of the optimization determine the tolerance of the convergence for the optimization. As a further research, various optimization methods can be performed while increasing number of kinematic parameters optimized searching different cost functions.

Notations

:Leg lengths
:Dir. along the x-axis of mobile platform
:Dir. along the y-axis of mobile platform
:Dir. along the z-axis of mobile platform
:Pitch angle of mobile platform
:Roll angle of mobile platform
:Yaw angle of mobile platform
:Condition number of Jacobian matrix
Jacobian matrix
:Mobile platform center of gravity
:Base platform center of gravity
:Angle between
:Angle between
:Radius of a circle passing through center of joints between base platform and legs
:Radius of a circle passing through center of joints between mobile platform and legs
:Angle between connection points of legs on base platform
:Angle between connection points of legs on mobile platform
Rotational transformation matrix.

Data Availability

The Matlab Codes data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that they have no conflicts of interest.