# Integrated Trajectory Planning and Sloshing Suppression for Three-Dimensional Motion of Liquid Container Transfer Robot Arm

Authors: Wisnu Aribowo; Takahito Yamashita; Kazuhiko Terashima
Journal: Journal of Robotics
 (2015)

---

## Abstract

For liquid transfer system in three-dimensional space, the use of multijoint robot arm provides much flexibility. To realize quick point-to-point motion with minimal sloshing in such system, we propose an integrated framework of trajectory planning and sloshing suppression. The robot motion is decomposed into translational motion of the robot wrist and rotational motion of the robot hand to ensure the upright orientation of the liquid container. The trajectory planning for the translational motion is based on cubic spline optimization with free via points that produces smooth trajectory in joint space while it still allows obstacle avoidance in task space. Input shaping technique is applied in the task space to suppress the motion induced sloshing, which is modeled as spherical pendulum with moving support. It has been found through simulations and experiments that the proposed approach is effective in generating quick motion with low amount of sloshing.


---

## Body

## 1. Introduction

In this research, we address a case where a robot arm is used to move a liquid filled container in three-dimensional space. Some example application cases for this setting are transfer of molten metal in casting industries and liquid containers handling by service robots. The use of articulated robot arm with high degree of freedom for these purposes provides more dexterity and flexibility when compared with the more dedicated type of liquid transfer devices. There are two main performance characteristics to consider in this system: quick motion and minimal sloshing. Those two criteria are conflicting, as quick motion tends to generate much sloshing of the contained liquid. Sloshing may cause the liquid to spill out of the container, and in the case of molten metal too much sloshing may affect the physical properties of the liquid. Moreover, other task, for example, pouring, cannot be performed when residual sloshing still occurs, thereby potentially reducing the system productivity. Thus the challenge in this problem is to generate trajectories with concurrent consideration of motion time minimization and sloshing suppression for point to point motion in a three-dimensional working space where static obstacles may exist.The motion speed of robot arm is mainly constrained by its kinematic and dynamic limit of the joints. One method to generate trajectories for such system is by posing the optimal control problem as a nonlinear optimization problem where the system states are discretized in some way. To reduce the large solution space and its associated computation time, the trajectory is often parameterized into a known structure, for example, polynomial function. Among the polynomial splines, cubic spline is the lowest order polynomial trajectory that guarantees continuation in acceleration and velocity as well as limited jerk. Higher order polynomials would result in smoother profile but longer motion time, while lower order ones would result in unlimited jerk.Numerous cubic spline algorithms and optimization schemes have been proposed in the literature [1–5]. Most of them try to optimize either the motion time or the location of knots, depending on the application needs. Relatively few addressed the optimization of time and location of knots concurrently [5]. In the latter approach, the location of knots is not fixed, which is suitable for trajectory planning problem with obstacle avoidance. In addition, it often results in a quicker motion because there is less restriction for the optimization algorithm to reach optimal velocity profile.The smooth cubic spline trajectories induce less vibration when compared to nonsmooth or lower order trajectories. However, in cases where minimal sloshing is required, smooth trajectory alone is not enough, and it needs to be supplemented with explicit sloshing control strategy. There have been numerous researches dealing with sloshing suppression in liquid container transfer system. Feedback control based approaches present good and relatively robust sloshing control ability [6, 7]. However, in some practical cases, sensor measurements cannot always be reliable or even are difficult to perform, for example, in the case of controlling sloshing of high-temperature molten metal. In these cases, feedforward approaches [8–14] are more useful. Based on sloshing models, the behavior of the system can be predicted, and thus the trajectory or motion path can be designed accordingly. Examples include generalized bang-bang control [9], infinite impulse response (IIR) filter [8], input shaping filter [10–13], and Hybrid Shape Approach filters [14]. Among them, the input shaping approach is particularly interesting because of its simplicity and effectiveness in practice. It works by decomposing the command input signal into two or more, using only the knowledge of natural frequency and damping ratio of the system [15, 16]. Over time, various type of input shapers with different characteristics have been proposed, and they have been successfully applied in practice for many kind of vibration problems.With respect to the motion path, traditionally sloshing is analyzed for one-dimensional straight motion, where the sloshing is often modeled as a simple pendulum. The motion planning and sloshing suppression is relatively simple and straightforward for this kind of problem [9, 13]. Researches that address suppression of sloshing or pendulum sway in higher dimensional space are relatively few. Tzamtzi et al. [17] modeled the sloshing as simple pendulum in 2D planar vertical motion. Williams et al. [18] used a dynamic programming approach to solve the sway-free motion of suspended spherical pendulum in 2D planar horizontal motion. Yano and Terashima [14] proposed HSA approach to control sloshing in three degree of freedom Cartesian liquid container transfer system.In this paper, we propose the integration of the trajectory planning of seven degrees of freedom liquid container transfer robot arm and the sloshing suppression of the contained liquid. We use the Mitsubishi PA 10-7C robot arm as the experiment device. The container is allowed to rotate about vertical axis while the pitch and roll angles are regulated to maintain upright orientation of the liquid container. Such system setup is different from that used in other related papers (e.g., [14, 17, 18]); in that translational motion in three-dimensional space plus orientation of the container needs to be controlled. The output of the proposed method is joint trajectories, to be used as a reference of the robot motion, thus obviating the need to change the closed-loop controller. Consequently, the proposed solution is easier to be applied to different kind of systems.Because the robot arm has spherical wrist, the motion can be conveniently decoupled into translational motion of the wrist locus and rotational motion of robot hand. This decomposition is useful in our application case, because for sloshing suppression we want to take full control of the container orientation. For the trajectory generation of the wrist translational motion, we define the problem as a nonlinear optimization model where each joint trajectory is parameterized as cubic spline. The locations of knots are part of decision variables, thus not fixed, to allow obstacle avoidance in task space and quick motion. It is an extension to a more basic cubic spline optimization reported in our previous paper [5] in two points: motion in three-dimensional space is addressed by using multijoint robot arm, and obstacle avoidance constraints are added. As for the sloshing suppression, we use a new strategy of applying input shaping. The decoupled translational and rotational motions are each shaped by using its respective suitable input shaper and then combined as the final trajectory. We use the spherical pendulum model to simulate the response of the trajectory. In addition, a numerical sloshing simulation model is built for the purpose of examining the sloshing behavior in response to vertical motion of the liquid container in 3D space. Finally, experiments are carried out for a few example cases.


## 2. Sloshing Model and Suppression

### 2.1. Equivalent Sloshing Model

The theoretical natural frequency of sloshing inside an upright cylindrical container can be derived from Navier-Stokes equations as [19](1)ωn2=gξnRtanh⁡ξnhR,where g is the gravity constant, R is the container radius, h is the liquid height, and ξn is the root of the derivative of Bessel function of the first kind where n is the sloshing mode. For the fundamental first mode, ξ1 equals 1.841.For control analysis, the phenomenon of sloshing is often modeled by its simpler mechanical equivalents, for example, pendulum and mass-spring-damper system. Here, we model the sloshing as a simple pendulum, where the sway of the pendulum corresponds to the elevation of the liquid surface. A simple pendulum system, as depicted in Figure1(a), consists of a point mass suspended by a massless rigid cable of length L to a friction-less movable support. For the case where the support moves in straight line, it can be described by the following equation of motion:(2)Lθ¨+gsin⁡θ=-x¨cos⁡⁡θ,where g is the gravity constant and x is the distance of the support from origin.Figure 1
(a) Simple pendulum and (b) spherical pendulum.
(a)
(b)In cases where the support moves in two-dimensional horizontal space, the system can be modeled as a spherical pendulum. Figure1(b) illustrates the spherical pendulum model, where the pendulum bob can freely swing in two degrees of freedom spherical space as response to movement of the support in horizontal planar space. Angle θ is the angle between planeZY and the mass, and φ is the angle in the perpendicular direction.The position of the pendulum bob in Cartesian coordinate system is(3)xp=-Lsin⁡θcos⁡⁡φ,yp=Lsin⁡θ,zp=-Lcos⁡⁡θcos⁡⁡φ.The equations of motion of the spherical pendulum can be derived as follows [18]:(4)θ¨=2θ˙φ˙tan⁡φ-gLsin⁡θcos⁡⁡φ+1Lcos⁡⁡θcos⁡⁡φx¨,φ¨=-θ˙2φ˙tan⁡φ-gLcos⁡⁡θsin⁡φ-1Lx¨sin⁡θsin⁡φ+y¨cos⁡⁡φ,where x and y are the position of the support. Indeed, when there is no acceleration inY direction, and initially φ=0, the above equations of motion will reduce to the simple pendulum case in (2).


### 2.2. Input Shaping

In the simple pendulum case, the natural frequency and damping ratio can be obtained by linearization of the equation of motion around equilibrium as(5)ω0=gL,ζ=c2mLg.When excited by an impulse, the sway response of the pendulum is decaying sinusoidal:(6)θt=Aω01-ζ2e-ζω0tsin⁡ω0t1-ζ2.Such sway can be suppressed by shaping the impulse input such that the sum of the response is zero. It is performed by splitting the original impulse into two or more impulses and delays them such that their responses cancel each other [15, 16]. For example, in the simplest case of splitting into two impulses, the amplitudes of the first and second impulses are K/1+K and 1/1+K, respectively, where (7)K=e-ζπ/1-ζ2,and time between impulse equals one-half damped period of vibration:(8)ΔT=πω01-ζ2.This set of impulses is called an input shaper, and the two-impulse input shaper above is known as zero vibration (ZV) shaper. A more robust input shaper type can be obtained by convolving two ZV input shapers together, resulting in three-impulse input shaper, known as zero vibration and derivative (ZVD) input shaper. Other than those basic input shapers, there are still many more input shapers developed with different characteristics and capabilities to tackle different type of vibration problems.In a linear system, any signal can be represented as a collection of scaled and shifted impulses. Therefore the system response to a signal is the linear sum of the responses of each component impulse. It follows that the input shaping can be applied by convolving the system input signal with the desired input shaper, as shown in Figure2. As a result, the motion time becomes longer, as much as the length of the input shaper. For example, applying ZV input shaper increases the motion time by half damped period of vibration, while with ZVD input shaper the motion time becomes one damped period longer. For general planar horizontal motion, the input signal may be in the form of acceleration or velocity reference of the motion.Figure 2
Shaping arbitrary signal by convolution.For systems whose motion is commanded by rotational commands, for example, rotary crane and other rotating mechanical systems, regular input shapers do not work well for suppressing motion induced vibration. Those problems exhibit different characteristics in that the angular command impulse excites sway in radial as well as tangential directions. As the angular velocity increases, the natural frequency bifurcates from the nominal natural frequency, resulting in two sway modes in the system. For this kind of problem, Lawrence and Singhose [20] developed ZV2lin shaper. The input shaper consists of three impulses, where the amplitude and timing are determined by solving an optimization problem which seeks minimization of vibration at the end of impulse train. Those are determined not only by natural frequency and damping ratio, but also by the nominal angular velocity and rotation radius. The input shaper length is approximately equal to one damped period, more or less comparable to the length of ZVD shaper.


## 2.1. Equivalent Sloshing Model

The theoretical natural frequency of sloshing inside an upright cylindrical container can be derived from Navier-Stokes equations as [19](1)ωn2=gξnRtanh⁡ξnhR,where g is the gravity constant, R is the container radius, h is the liquid height, and ξn is the root of the derivative of Bessel function of the first kind where n is the sloshing mode. For the fundamental first mode, ξ1 equals 1.841.For control analysis, the phenomenon of sloshing is often modeled by its simpler mechanical equivalents, for example, pendulum and mass-spring-damper system. Here, we model the sloshing as a simple pendulum, where the sway of the pendulum corresponds to the elevation of the liquid surface. A simple pendulum system, as depicted in Figure1(a), consists of a point mass suspended by a massless rigid cable of length L to a friction-less movable support. For the case where the support moves in straight line, it can be described by the following equation of motion:(2)Lθ¨+gsin⁡θ=-x¨cos⁡⁡θ,where g is the gravity constant and x is the distance of the support from origin.Figure 1
(a) Simple pendulum and (b) spherical pendulum.
(a)
(b)In cases where the support moves in two-dimensional horizontal space, the system can be modeled as a spherical pendulum. Figure1(b) illustrates the spherical pendulum model, where the pendulum bob can freely swing in two degrees of freedom spherical space as response to movement of the support in horizontal planar space. Angle θ is the angle between planeZY and the mass, and φ is the angle in the perpendicular direction.The position of the pendulum bob in Cartesian coordinate system is(3)xp=-Lsin⁡θcos⁡⁡φ,yp=Lsin⁡θ,zp=-Lcos⁡⁡θcos⁡⁡φ.The equations of motion of the spherical pendulum can be derived as follows [18]:(4)θ¨=2θ˙φ˙tan⁡φ-gLsin⁡θcos⁡⁡φ+1Lcos⁡⁡θcos⁡⁡φx¨,φ¨=-θ˙2φ˙tan⁡φ-gLcos⁡⁡θsin⁡φ-1Lx¨sin⁡θsin⁡φ+y¨cos⁡⁡φ,where x and y are the position of the support. Indeed, when there is no acceleration inY direction, and initially φ=0, the above equations of motion will reduce to the simple pendulum case in (2).


## 2.2. Input Shaping

In the simple pendulum case, the natural frequency and damping ratio can be obtained by linearization of the equation of motion around equilibrium as(5)ω0=gL,ζ=c2mLg.When excited by an impulse, the sway response of the pendulum is decaying sinusoidal:(6)θt=Aω01-ζ2e-ζω0tsin⁡ω0t1-ζ2.Such sway can be suppressed by shaping the impulse input such that the sum of the response is zero. It is performed by splitting the original impulse into two or more impulses and delays them such that their responses cancel each other [15, 16]. For example, in the simplest case of splitting into two impulses, the amplitudes of the first and second impulses are K/1+K and 1/1+K, respectively, where (7)K=e-ζπ/1-ζ2,and time between impulse equals one-half damped period of vibration:(8)ΔT=πω01-ζ2.This set of impulses is called an input shaper, and the two-impulse input shaper above is known as zero vibration (ZV) shaper. A more robust input shaper type can be obtained by convolving two ZV input shapers together, resulting in three-impulse input shaper, known as zero vibration and derivative (ZVD) input shaper. Other than those basic input shapers, there are still many more input shapers developed with different characteristics and capabilities to tackle different type of vibration problems.In a linear system, any signal can be represented as a collection of scaled and shifted impulses. Therefore the system response to a signal is the linear sum of the responses of each component impulse. It follows that the input shaping can be applied by convolving the system input signal with the desired input shaper, as shown in Figure2. As a result, the motion time becomes longer, as much as the length of the input shaper. For example, applying ZV input shaper increases the motion time by half damped period of vibration, while with ZVD input shaper the motion time becomes one damped period longer. For general planar horizontal motion, the input signal may be in the form of acceleration or velocity reference of the motion.Figure 2
Shaping arbitrary signal by convolution.For systems whose motion is commanded by rotational commands, for example, rotary crane and other rotating mechanical systems, regular input shapers do not work well for suppressing motion induced vibration. Those problems exhibit different characteristics in that the angular command impulse excites sway in radial as well as tangential directions. As the angular velocity increases, the natural frequency bifurcates from the nominal natural frequency, resulting in two sway modes in the system. For this kind of problem, Lawrence and Singhose [20] developed ZV2lin shaper. The input shaper consists of three impulses, where the amplitude and timing are determined by solving an optimization problem which seeks minimization of vibration at the end of impulse train. Those are determined not only by natural frequency and damping ratio, but also by the nominal angular velocity and rotation radius. The input shaper length is approximately equal to one damped period, more or less comparable to the length of ZVD shaper.


## 3. Trajectory Planning

The liquid transfer system concerned uses the Mitsubishi PA 10-7C robot arm. The frames placement and link dimension of the robot arm are shown in Figure3. Table 1 lists the limit value of each joint angle and angular velocity. A cylindrical liquid container is attached to an extension at the robot tip. The diameter and height of the liquid container are 0.15 mm and 0.25 m, respectively. The distance between robot tip and the center of the liquid container is 0.12 meter; thus the length of the last link is practically 0.2 meter. For our current application, one of the arm joints, the S3 joint, is locked so that only six joints are actively used. Table 2 lists the Denavit-Hartenberg parameters of the robot arm configuration. The transformation matrix of a frame relative to the preceding frame is as follows:(9)Ai=​i-1Ri​i-1pi01=cos⁡⁡θi-sin⁡θicos⁡⁡αisin⁡θisin⁡αiaicos⁡⁡θisin⁡θicos⁡⁡θicos⁡⁡αi-cos⁡⁡θisin⁡αiaisin⁡θi0sin⁡αicos⁡⁡αidi0001.Table 1
Limit values of joint velocities and angles.
Limit
Joint
S1
S2
S3
E1
E2
W1
W2
V
MAX [rad/s]
1
1
2
2
2
π
2
π
2
π
q
UB [degree]
177
91
174
137
255
165
360Table 2
Denavit-Hartenberg parameters.
Parameter
Joint
S1
S2
S3
E1
E2
W1
W2
α [rad]
π/2
−π/2
π/2
−π/2
π/2
−π/2
0
d [m]
0.315
0
0.45
0
0.5
0
0.2
a [m]
0
0
0
0
0
0
0Figure 3
(a) Frames placement and (b) dimension of the liquid container system and the base frame.
(a)
Frames placement
(b)
Dimension of links and liquid container (in meter)The position and orientation of a link with respect to other link can be obtained by chain product of the transformation matrices:(10)​iTj=Ai+1Ai+2⋯Aj=Rij​ipj01.The upper left3×3 rotation matrix Rij describes the link orientation, while the vector pij represents the position.The links E2, W1, and W2 form a spherical wrist, where the three rotational joint axes intersect at a common point. This robot configuration simplifies the inverse kinematic calculation because it can be decoupled into two simpler subsystems. The first three rotational links S1, S2, and E1 provide the positioning in three-dimensional space, while the spherical wrist adjusts the orientation of end effector.In the inverse kinematics for position, the position of wrist locus in the task space (i.e.xw, yw, and zw) is known, and the associated configuration of S1, S2, and E1 joints (i.e. θ1, θ2, and θ4) are to be sought. The values can be calculated geometrically as follows:(11)θ1=atan⁡2yw,xw,θ4=atan⁡21-D2,D,θ2=atan⁡2zw-L1,xw2+yw2-atan⁡2L3sin⁡θ4,L2+L3cos⁡⁡θ4-μ2,where Li is the length of link i, atan2 is the arc tangent function with two arguments, and(12)D=xw2+yw2+zw-L12-L22-L322L2L3.Meanwhile, the inverse kinematics for orientation calculates the configuration of E2, W1, and W2 joints (i.e.θ5, θ6, and θ7) when the desired orientation of end effector relative to the base frame R70 is known. First, from the inverse kinematics for position, the rotation matrix describing the orientation of the third link can be obtained R40. Then, the rotation operation of the wrist joints required to realize the end effector posture is calculated as(13)R74=R40TR70.The three Euler angles can be calculated from the matrix, and finally we can obtain the configuration of the wrist joints:(14)θ5=atan⁡2R742,3,R741,3±π,(15)θ6=atan⁡21-R473,32,R743,3,(16)θ7=atan⁡2R743,2,R743,1±π,where R74r,c is the element of the matrix R74 at row r and column c. When there are more than one possibilities of angle values, choose the one suitable for the application.Accordingly, the trajectory planning of the liquid container transfer system is decoupled into planning of the translation of the robot wrist position and rotation of the robot hand. The translation part is designed as joint space cubic spline trajectory, while the hand rotation is designed in task space to maintain upright orientation of the liquid container. Upright orientation means that the container has only one degree of freedom: the rotation along the verticalz-axis of fixed coordinate system. Those translation and rotation trajectories are then combined to make the final quick slosh-free trajectory solution. The trajectory planning steps are shown in Figure 4, which clearly depicts whether each step is done in joint space or task space. Transformation between joint and task space is carried out using direct and inverse kinematics.Figure 4
Trajectory planning steps for (a) translational motion and (b) hand orientation.
(a)
Trajectory planning for translation motion
(b)
Trajectory planning for hand orientationFigure4(a) shows the overview of planning the translational motion. It starts with specifying the start point and end point, as well as several knots between them in task space, which are then transformed to joint space as the initial solution to the cubic spline optimization. The optimization step encompasses joint space as well as task space. Although the trajectory is generated in joint space, it is possible to define task space constraints in the optimization, for example, obstacle avoidance. The obtained solution is joint space trajectory, which is transformed back to task space, where the command shaping takes place. The solution of the translational motion is obtained as joint trajectories of links S1, S2, and E1.Figure4(b) shows the steps of planning the rotational motion of the robot hand. The input is the start and end angle of the container alongz-axis relative toXZ vertical plane. The shaped rotation trajectory is then transformed to joint space as the solution of the rotation part of the motion in form of joint trajectories of links E2, W1, and W2. Although rotation trajectory generation is independent of translation trajectory generation, it has to be performed after because the former needs the orientation history of the translational motion as well as the total motion time in order to coordinate the overall motion.


### 3.1. Translational Motion Planning

In this section, we outline the cubic spline optimization approach for trajectory planning. It builds from the same basic approach as in a previous paper [5] for two-dimensional problem and then adapted to the current application of motion in three-dimensional space and is complemented with input shaping in the task space.A cubic spline trajectory of a robot joint is comprised of several curve segments, as illustrated in Figure5, where each segment is represented by a cubic polynomial function of time. The general form of the cubic function in a segment i is(17)qit=ait-ti3+bit-ti2+cit-ti+di,where t denotes the time, which ranges from ti to ti+1.Figure 5
Cubic spline trajectory.The curve segments are connected by knots. Thus for a trajectory that consists ofn curve segments, there are n+1 prespecified joint angle values θi, i=1⋯n+1 which consist of a start point, an end point, and n-1 knots between them. The trajectory is required to pass through all those points. Other requirements are continuity in velocity and acceleration at every knot. The relationship between joint angle θi, joint velocity θ˙i, joint acceleration θ¨i, and time interval hi=ti+1-ti is as follows: (18)θ¨ihi+2θ¨i+1hi+hi+1+θ¨i+2hi+1=6θi+2-θi+1hi+1-θi+1-θihi.In a joint trajectory withn segments (or n+1 points), there are n-1 such equations. Together with the boundary requirements of acceleration (θ¨1=θ¨n+1=0), we have a total of n+1 equations. The robot motion usually also requires that velocities at start and end points equal zero. To accommodate this, two more unknowns are added by inserting two virtual points: one after the start point and another one just before the end point. Furthermore, because θ¨1 and θ¨n+1 are already known as zero, they can be removed from the equations. Joint angles θi are defined beforehand by the trajectory designer. The knots may be specified as points in Cartesian task space and later transformed into joint space by inverse kinematics. The above system of linear equations forms a symmetric tridiagonal system, which can be solved efficiently using the tridiagonal matrix algorithm (TDMA) to find acceleration at each knot.The trajectories of S1, S2, and E1 joints are each designed as piecewise cubic splines. The joint trajectories are coupled to each other by using the same value of time interval variableshi. We formulate the trajectory generation problem as a nonlinear optimization problem. The objective function is combination of motion time and weighted squared acceleration. There are several constraints to consider in the trajectory generation: robot position in the task space, joint velocity, and joint jerk. According to (18), there are three kinds of variables involved: segment time, joint acceleration, and joint angle at knots. We let both the segment time and the joint angles (knots location) as the optimization variables. The optimization problem is thus as follows:(19)min⁡xF=∑i=1nhi+wa∫q¨i2dtsubjecttoqij∗≤qjMAX;i=1⋯n,j=1⋯mvij∗≤vjMAX;i=1⋯n,j=1⋯mθ⃛≤JMAX;i=1⋯n+1,j=1⋯mdiLj,O≥rLj+rO;i=1⋯n,j=1⋯m,where x is the vector of optimization variables which contains the segment times as well as the joint angles excluding the fixed start point, end point, and the two virtual points:(20)x=h1,…,hn,θij,…,θnm;i=3⋯n-1,j=1⋯m.Here, all internal knots are part of the optimization variables instead of predefined points. This provides flexibility for the optimization to generate quick motion as well as obstacle avoidance at the same time. If the planned motion is required to pass certain fixed points on the way, they can be added as additional equality constraints. The first term of the objective function is to minimize the total motion time, which is the main objective. The second term is used to avoid unnecessary motion of robot joints, wherewa is the weight parameter. In our subsequent simulations and experiments, the parameter value equals 1e-4.The first and second constraints are joint angle and velocity constraints, according to the physical joint limits in Table1. The third is jerk constraint, to prevent excessive value of joint jerk. The maximum jerk of each joint is set equal to 15 m/s3. The last constraint is for obstacle avoidance, where the obstacle and robot links are modeled as rigid capsules with radius rO and rLj, respectively. The function diLj,O is the minimum distance between link Lj and static obstacle O in segment i. The obstacle and links’ axis are discretized into several representative points. The distance between a link and obstacle is defined as the minimum distance between the points in the obstacle and the points in the link.We use the MATLAB Optimization Toolbox implementation of Sequential Quadratic Programming (SQP) for solving the above constrained nonlinear trajectory planning optimization problem. The SQP is one of the most used methods in solving general constrained nonlinear optimization problems. The method solves the problem by iteratively solving a quadratic subproblem, which is the quadratic approximation of the Lagrange function of the original nonlinear problem. The solution of the quadratic subproblem is then used to form a new set of optimization variables for the next iteration of the original problem. The calculation time differs depending on the computer used for calculation; it is around 15 seconds for a typical problem size in our Core2 Duo computer.Figure6 shows an example motion path of a cubic trajectory generated by the above method in perspective view and top view. The PA-10 robot arm is stationed at origin of the world frame. The robot has to move its wrist from start point [0.3,-0.8,0.2] to end point [0.5,0.5,0.2] in the task space, which correspond to θSTART=-1.21,-1.25,-0.87 rad and θEND=0.79,-0.97,-1.44 rad in the joint space. The robot posture shown in the figure is the configuration at start point. A cylindrical obstacle with radius 0.05 m lies vertically in the workspace at [0.5,-0.3]. For this problem, four initial via points are specified randomly between the start and end points. The line going from the start point and the end point is the generated motion path in the task space, which successfully avoids the obstacle. The small spheres on the line denote the location of output via points. The total motion time is 2.7 seconds. Figure 7 shows the joint trajectories of the example case: the piecewise cubic joint angles and the piecewise quadratic joint velocities. The small circles in the figure denote the respective values at the knots, including at the two virtual knots. From the velocity figure, we can understand that the motion time is limited by joint S1, where the velocity hits its limit most of the time.Figure 6
Motion path of example case of translation motion.
(a)
Perspective view
(b)
Top view (X-Y plane)Figure 7
Trajectory of example case of translational motion: (a) joint angle (b) joint velocity.
(a)
Joint angle
(b)
Joint velocityThe next step in the translational motion planning is to suppress the sloshing generated by the trajectory obtained above. For analysis, we use the pendulum model in Section2.1 to represent the sloshing, where the pendulum support is attached at the robot wrist locus and the pendulum sways freely as the robot moves. For the case where the pendulum support can move in a horizontal planar space, we first decompose the input signal to its projection into orthogonal Cartesian axes and then apply input shaping independently to each input signal. By doing this, the spherical pendulum is considered as two independent simple pendulums which act on orthogonal vertical planes such that the sway behavior in one subsystem is not affected by the states in the other subsystem. This simplified approach works well in small angle regions, and indeed linearization of (4) by small angle approximations and removing quadratic and higher terms results in two independent equations of simple pendulum.We use the ZVD input shaping for sloshing suppression. For that, the velocity trajectory is decomposed into two orthogonal velocities parallel tox-axis andy-axis, respectively. Each velocity trajectory is then convolved individually with ZVD input shaper. The shaped results are then added back and integrated to produce the shaped trajectory in task space. Joint space solution can be obtained by inverse kinematics.As an example, we simulate a spherical pendulum of length 0.2 meter moving according to the trajectory obtained above. The natural frequency and the vibration period of the pendulum model are 7.0 rad/sec and 0.8973 seconds, respectively. In sloshing problems, the damping is usually very small; therefore here we assume zero damping in the pendulum model. Figure8 shows the original and ZVD shaped Cartesian velocities inX andY directions. The motion time is 2.7 seconds and 3.6 seconds in the unshaped and ZVD shaped case, respectively. Figure 9 shows the sway angles in the spherical coordinate system. At the end of motion, sway is much suppressed in the shaped trajectory. Figure 10 compares the projection of the pendulum position at horizontal plane for the two cases: original and shaped trajectories. It is clear from the figure that the shaped trajectory generates significantly less sway than its unshaped counterpart.Figure 8
Original and ZVD shaped velocity in task space: (a)X direction (b)Y direction.
(a)
X direction
(b)
Y directionFigure 9
Pendulum sway angles: (a) without shaping (b) with ZVD shaping.
(a)
Without shaping
(b)
ZVD input shapingFigure 10
Projection of pendulum bob position inXY plane: (a) without shaping (b) with ZVD shaping.
(a)
Without shaping
(b)
ZVD input shaping


### 3.2. Rotational Motion Planning

The liquid container has to be maintained upright all time. This means that there is only one task space DOF: the rotation along the verticalz-axis. Thus the hand rotational motion can be simply parameterized by the angle alongz-axis relative to theXZ plane. Here, robot hand means any parts after the wrist, including the last robot link, the container holder, and the liquid container itself. The length of the robot hand is the distance between the wrist locus and the center of container, which equals 0.2 meter. The trajectory is generated by applying ZV2lin input shaping to a rectangular velocity profile, in which the maximum angular velocity is (21)s˙=γf-γsT,where T is the time of the unshaped translation motion and γs and γf are the desired start and finish angle value, respectively. This is to coordinate the translational and rotational motions, so that both motions end at approximately the same time.As an example, we use the same spherical pendulum model to simulate the motion of the pendulum support on a circular trajectory of lengthπ/2 radian, where the robot hand moves from angle -π/2 radian to 0 (relative toXZ plane). The unshaped translation motion time is 2.7 seconds; thus the maximum angular velocity for the rotational motion equals 0.58 rad/sec. Upon applying ZV2lin input shaper, the rotational motion time becomes 3.47 seconds.Figure11 shows the comparison of the projections of the pendulum bob onXY plane of the unshaped trajectory and the trajectory shaped by ZV2lin input shaper. The circular motion starts from the origin and stops at (0.2, 0.2). In the shaped case, the pendulum sway is much suppressed, both the sway while the support is moving and the residual sway after the support stops. This shows the effectiveness of the shaper in suppressing vibration for motions generated by angular command.Figure 11
Pendulum position inXY plane of circular motion for (a) without shaping (b) with shaping.
(a)
Without shaping
(b)
ZV2lin shapedThe trajectories of the three spherical wrist links can be obtained using the inverse kinematic relations in (13)–(16), the shaped angle value, and the orientation history of the translation motion part. The rotation matrix R70 is the desired orientation of the liquid container:(22)R70=cos⁡⁡γ-sin⁡γ0sin⁡γcos⁡⁡γ0001,where γ is the angle aboutz-axis.The translational and rotational motions are independent as each is handled by different set of robot joints. Because of that, the final trajectory can be obtained simply by combining the joint trajectories of the translational and rotational motions. Once again, we simulate cases where spherical pendulum is attached to the robot tip, where the robot wrist locus moves from location (0.3, −0.8, 0.2) to (0.5, 0.5, 0.2) while at the same time the robot hand rotates with respect to the wrist from angle-π/2 rad to 0 rad. In short, the trajectory is the combination of the simulated trajectories at Sections 3.1 and 3.2. Figure 12 compares the sway angles of the case without shaping and the shaped case. The shaping is done individually: the translation motion is shaped using ZVD input shaper and the rotational motion is shaped using ZV2lin input shaping. The motion time is the same as in Section 3.1: 2.7 seconds and 3.6 seconds for the unshaped and shaped case, respectively. Figure 13 shows the motion path of the robot tip and the projection of the pendulum bob inXY plane. From those two figures, we can understand that the strategy of combining individually shaped trajectories result in much reduced sway of the pendulum. In that way, for different conditions (in this case, translational and rotational motion), we can use suitable input shaping for each condition and then combined the generated trajectories.Figure 12
Pendulum sway angles: (a) without shaping (b) combined input shaping.
(a)
Without shaping
(b)
Combined input shapingFigure 13
Projection of pendulum bob position inXY plane: (a) without shaping (b) combined input shaping.
(a)
Without shaping
(b)
Combined input shaping


## 3.1. Translational Motion Planning

In this section, we outline the cubic spline optimization approach for trajectory planning. It builds from the same basic approach as in a previous paper [5] for two-dimensional problem and then adapted to the current application of motion in three-dimensional space and is complemented with input shaping in the task space.A cubic spline trajectory of a robot joint is comprised of several curve segments, as illustrated in Figure5, where each segment is represented by a cubic polynomial function of time. The general form of the cubic function in a segment i is(17)qit=ait-ti3+bit-ti2+cit-ti+di,where t denotes the time, which ranges from ti to ti+1.Figure 5
Cubic spline trajectory.The curve segments are connected by knots. Thus for a trajectory that consists ofn curve segments, there are n+1 prespecified joint angle values θi, i=1⋯n+1 which consist of a start point, an end point, and n-1 knots between them. The trajectory is required to pass through all those points. Other requirements are continuity in velocity and acceleration at every knot. The relationship between joint angle θi, joint velocity θ˙i, joint acceleration θ¨i, and time interval hi=ti+1-ti is as follows: (18)θ¨ihi+2θ¨i+1hi+hi+1+θ¨i+2hi+1=6θi+2-θi+1hi+1-θi+1-θihi.In a joint trajectory withn segments (or n+1 points), there are n-1 such equations. Together with the boundary requirements of acceleration (θ¨1=θ¨n+1=0), we have a total of n+1 equations. The robot motion usually also requires that velocities at start and end points equal zero. To accommodate this, two more unknowns are added by inserting two virtual points: one after the start point and another one just before the end point. Furthermore, because θ¨1 and θ¨n+1 are already known as zero, they can be removed from the equations. Joint angles θi are defined beforehand by the trajectory designer. The knots may be specified as points in Cartesian task space and later transformed into joint space by inverse kinematics. The above system of linear equations forms a symmetric tridiagonal system, which can be solved efficiently using the tridiagonal matrix algorithm (TDMA) to find acceleration at each knot.The trajectories of S1, S2, and E1 joints are each designed as piecewise cubic splines. The joint trajectories are coupled to each other by using the same value of time interval variableshi. We formulate the trajectory generation problem as a nonlinear optimization problem. The objective function is combination of motion time and weighted squared acceleration. There are several constraints to consider in the trajectory generation: robot position in the task space, joint velocity, and joint jerk. According to (18), there are three kinds of variables involved: segment time, joint acceleration, and joint angle at knots. We let both the segment time and the joint angles (knots location) as the optimization variables. The optimization problem is thus as follows:(19)min⁡xF=∑i=1nhi+wa∫q¨i2dtsubjecttoqij∗≤qjMAX;i=1⋯n,j=1⋯mvij∗≤vjMAX;i=1⋯n,j=1⋯mθ⃛≤JMAX;i=1⋯n+1,j=1⋯mdiLj,O≥rLj+rO;i=1⋯n,j=1⋯m,where x is the vector of optimization variables which contains the segment times as well as the joint angles excluding the fixed start point, end point, and the two virtual points:(20)x=h1,…,hn,θij,…,θnm;i=3⋯n-1,j=1⋯m.Here, all internal knots are part of the optimization variables instead of predefined points. This provides flexibility for the optimization to generate quick motion as well as obstacle avoidance at the same time. If the planned motion is required to pass certain fixed points on the way, they can be added as additional equality constraints. The first term of the objective function is to minimize the total motion time, which is the main objective. The second term is used to avoid unnecessary motion of robot joints, wherewa is the weight parameter. In our subsequent simulations and experiments, the parameter value equals 1e-4.The first and second constraints are joint angle and velocity constraints, according to the physical joint limits in Table1. The third is jerk constraint, to prevent excessive value of joint jerk. The maximum jerk of each joint is set equal to 15 m/s3. The last constraint is for obstacle avoidance, where the obstacle and robot links are modeled as rigid capsules with radius rO and rLj, respectively. The function diLj,O is the minimum distance between link Lj and static obstacle O in segment i. The obstacle and links’ axis are discretized into several representative points. The distance between a link and obstacle is defined as the minimum distance between the points in the obstacle and the points in the link.We use the MATLAB Optimization Toolbox implementation of Sequential Quadratic Programming (SQP) for solving the above constrained nonlinear trajectory planning optimization problem. The SQP is one of the most used methods in solving general constrained nonlinear optimization problems. The method solves the problem by iteratively solving a quadratic subproblem, which is the quadratic approximation of the Lagrange function of the original nonlinear problem. The solution of the quadratic subproblem is then used to form a new set of optimization variables for the next iteration of the original problem. The calculation time differs depending on the computer used for calculation; it is around 15 seconds for a typical problem size in our Core2 Duo computer.Figure6 shows an example motion path of a cubic trajectory generated by the above method in perspective view and top view. The PA-10 robot arm is stationed at origin of the world frame. The robot has to move its wrist from start point [0.3,-0.8,0.2] to end point [0.5,0.5,0.2] in the task space, which correspond to θSTART=-1.21,-1.25,-0.87 rad and θEND=0.79,-0.97,-1.44 rad in the joint space. The robot posture shown in the figure is the configuration at start point. A cylindrical obstacle with radius 0.05 m lies vertically in the workspace at [0.5,-0.3]. For this problem, four initial via points are specified randomly between the start and end points. The line going from the start point and the end point is the generated motion path in the task space, which successfully avoids the obstacle. The small spheres on the line denote the location of output via points. The total motion time is 2.7 seconds. Figure 7 shows the joint trajectories of the example case: the piecewise cubic joint angles and the piecewise quadratic joint velocities. The small circles in the figure denote the respective values at the knots, including at the two virtual knots. From the velocity figure, we can understand that the motion time is limited by joint S1, where the velocity hits its limit most of the time.Figure 6
Motion path of example case of translation motion.
(a)
Perspective view
(b)
Top view (X-Y plane)Figure 7
Trajectory of example case of translational motion: (a) joint angle (b) joint velocity.
(a)
Joint angle
(b)
Joint velocityThe next step in the translational motion planning is to suppress the sloshing generated by the trajectory obtained above. For analysis, we use the pendulum model in Section2.1 to represent the sloshing, where the pendulum support is attached at the robot wrist locus and the pendulum sways freely as the robot moves. For the case where the pendulum support can move in a horizontal planar space, we first decompose the input signal to its projection into orthogonal Cartesian axes and then apply input shaping independently to each input signal. By doing this, the spherical pendulum is considered as two independent simple pendulums which act on orthogonal vertical planes such that the sway behavior in one subsystem is not affected by the states in the other subsystem. This simplified approach works well in small angle regions, and indeed linearization of (4) by small angle approximations and removing quadratic and higher terms results in two independent equations of simple pendulum.We use the ZVD input shaping for sloshing suppression. For that, the velocity trajectory is decomposed into two orthogonal velocities parallel tox-axis andy-axis, respectively. Each velocity trajectory is then convolved individually with ZVD input shaper. The shaped results are then added back and integrated to produce the shaped trajectory in task space. Joint space solution can be obtained by inverse kinematics.As an example, we simulate a spherical pendulum of length 0.2 meter moving according to the trajectory obtained above. The natural frequency and the vibration period of the pendulum model are 7.0 rad/sec and 0.8973 seconds, respectively. In sloshing problems, the damping is usually very small; therefore here we assume zero damping in the pendulum model. Figure8 shows the original and ZVD shaped Cartesian velocities inX andY directions. The motion time is 2.7 seconds and 3.6 seconds in the unshaped and ZVD shaped case, respectively. Figure 9 shows the sway angles in the spherical coordinate system. At the end of motion, sway is much suppressed in the shaped trajectory. Figure 10 compares the projection of the pendulum position at horizontal plane for the two cases: original and shaped trajectories. It is clear from the figure that the shaped trajectory generates significantly less sway than its unshaped counterpart.Figure 8
Original and ZVD shaped velocity in task space: (a)X direction (b)Y direction.
(a)
X direction
(b)
Y directionFigure 9
Pendulum sway angles: (a) without shaping (b) with ZVD shaping.
(a)
Without shaping
(b)
ZVD input shapingFigure 10
Projection of pendulum bob position inXY plane: (a) without shaping (b) with ZVD shaping.
(a)
Without shaping
(b)
ZVD input shaping


## 3.2. Rotational Motion Planning

The liquid container has to be maintained upright all time. This means that there is only one task space DOF: the rotation along the verticalz-axis. Thus the hand rotational motion can be simply parameterized by the angle alongz-axis relative to theXZ plane. Here, robot hand means any parts after the wrist, including the last robot link, the container holder, and the liquid container itself. The length of the robot hand is the distance between the wrist locus and the center of container, which equals 0.2 meter. The trajectory is generated by applying ZV2lin input shaping to a rectangular velocity profile, in which the maximum angular velocity is (21)s˙=γf-γsT,where T is the time of the unshaped translation motion and γs and γf are the desired start and finish angle value, respectively. This is to coordinate the translational and rotational motions, so that both motions end at approximately the same time.As an example, we use the same spherical pendulum model to simulate the motion of the pendulum support on a circular trajectory of lengthπ/2 radian, where the robot hand moves from angle -π/2 radian to 0 (relative toXZ plane). The unshaped translation motion time is 2.7 seconds; thus the maximum angular velocity for the rotational motion equals 0.58 rad/sec. Upon applying ZV2lin input shaper, the rotational motion time becomes 3.47 seconds.Figure11 shows the comparison of the projections of the pendulum bob onXY plane of the unshaped trajectory and the trajectory shaped by ZV2lin input shaper. The circular motion starts from the origin and stops at (0.2, 0.2). In the shaped case, the pendulum sway is much suppressed, both the sway while the support is moving and the residual sway after the support stops. This shows the effectiveness of the shaper in suppressing vibration for motions generated by angular command.Figure 11
Pendulum position inXY plane of circular motion for (a) without shaping (b) with shaping.
(a)
Without shaping
(b)
ZV2lin shapedThe trajectories of the three spherical wrist links can be obtained using the inverse kinematic relations in (13)–(16), the shaped angle value, and the orientation history of the translation motion part. The rotation matrix R70 is the desired orientation of the liquid container:(22)R70=cos⁡⁡γ-sin⁡γ0sin⁡γcos⁡⁡γ0001,where γ is the angle aboutz-axis.The translational and rotational motions are independent as each is handled by different set of robot joints. Because of that, the final trajectory can be obtained simply by combining the joint trajectories of the translational and rotational motions. Once again, we simulate cases where spherical pendulum is attached to the robot tip, where the robot wrist locus moves from location (0.3, −0.8, 0.2) to (0.5, 0.5, 0.2) while at the same time the robot hand rotates with respect to the wrist from angle-π/2 rad to 0 rad. In short, the trajectory is the combination of the simulated trajectories at Sections 3.1 and 3.2. Figure 12 compares the sway angles of the case without shaping and the shaped case. The shaping is done individually: the translation motion is shaped using ZVD input shaper and the rotational motion is shaped using ZV2lin input shaping. The motion time is the same as in Section 3.1: 2.7 seconds and 3.6 seconds for the unshaped and shaped case, respectively. Figure 13 shows the motion path of the robot tip and the projection of the pendulum bob inXY plane. From those two figures, we can understand that the strategy of combining individually shaped trajectories result in much reduced sway of the pendulum. In that way, for different conditions (in this case, translational and rotational motion), we can use suitable input shaping for each condition and then combined the generated trajectories.Figure 12
Pendulum sway angles: (a) without shaping (b) combined input shaping.
(a)
Without shaping
(b)
Combined input shapingFigure 13
Projection of pendulum bob position inXY plane: (a) without shaping (b) combined input shaping.
(a)
Without shaping
(b)
Combined input shaping


## 4. Numerical Simulation of Sloshing

In order to analyze the effect of vertical motion to sloshing, we built a numerical simulation of sloshing inside a moving liquid container based on exact distributed parameter computational fluid dynamics (CFD) model, not the simple pendulum model. The liquid container is cylindrical, with 150 mm diameter and 250 mm in height. The depth of liquid is 170 mm. The simulation model is developed in the open source OpenFOAM software package. The liquid oscillation is measured during simulation time with measurement sampling time 0.01 seconds by two virtual probes that record the liquid height near the container wall at the direction ofx-axis andy-axis, respectively. By using the measured sloshing data, we identified the dominant fundamental mode of the sloshing at frequency 2.3 Hz.As explained in Section3.1, input shaping for translational motion is implemented by decomposing the trajectory into two Cartesian components (x-axis andy-axis) and then applying input shaping to those two trajectory components individually. The shaped trajectories are then combined as the final shaped trajectory. By using the spherical pendulum model, it has been shown that the approach works well in suppressing the pendulum sway. However, when we consider motion in 3D space, the motion component inz-axis also has to be taken into consideration. To examine the effect of acceleration inz-axis toward sloshing, we run a simulation comparing trajectories with and without vertical acceleration. In both trajectories, the liquid container moves horizontally for the first one second to generate some sloshing. Then, in trajectory A, we let the container stay still, while in trajectory B, the liquid container moves upward. Figure 14 shows the comparison of the sloshing generated by both trajectories. From second one, the sloshing frequency differs slightly, where it is a bit higher when the liquid container accelerates upward. Figure 15 shows the comparison in frequency domain. This simulation result agrees with the theory that states natural frequency of sloshing is directly proportional to square root of vertical acceleration (1). In the usual case of lateral movement, it consists of only the constant gravitational acceleration. But when it accelerates upward, the vertical acceleration would be higher than g, and as a result natural frequency becomes higher. The theoretical natural frequency under constant vertical acceleration of 1 m/s2 equals 2.4 Hz, which is approximately equal with the simulation result.Figure 14
Comparison of the simulated sloshing under vertical acceleration.Figure 15
Comparison of the sloshing frequency under vertical acceleration.The consequence of this to the sloshing suppression is that we have to consider the frequency shift in the input shaper calculation. From simulation, we understand that the range of the frequency shift is not wide. With current setup of the robot arm, the vertical acceleration of a typical motion is under ±4 m/s2, which approximately corresponds to maximum frequency shift of ±0.4 Hz. To deal with this small shift, we use zero vibration and damping (ZVD) input shaper, which is a more robust version of the input shaper, instead of the basic ZV shaper.To examine the sloshing in a vertical dominant motion trajectory, we setup another simulation case where liquid container is moved from location (0.2, 0.6, 0.2) to location (0.3, 0.2, 1.1). In this trajectory, vertical movement is relatively dominant. Thus, we expect the natural frequency to shift along the trajectory. Figure16 shows the sloshing comparison. Consistent with the previous example case, the shaped trajectories result in less vibration. In addition, the robust ZVD shaped trajectory improves the sloshing suppression. The maximum residual sloshing is 46.22 mm, 4.04 mm, and 1.76 mm for original, ZV shaped, and ZVD shaped trajectories, respectively. From this simulation study, we understand that vertical motion affects the sloshing in that small shift of natural frequency occurs. The small shift justifies the use of the zero vibration and derivative (ZVD) input shaping, which is more robust than the default ZV input shaping, in our proposed solution.Figure 16
Comparison of sloshing in vertically dominant case.


## 5. Experiment Results and Discussion

In order to assess the effectiveness of the overall proposed system, including the trajectory planning for translational and rotational motion as well as the sloshing, we setup a few experiment cases. Table3 lists the start and end position and orientation. The start points and end points are the position of the center of the container liquid. The angles are measured relative toXZ plane. Figure 17 illustrates the motion path of the three experiment cases, where the start and end points of each case are shown by small spheres.Table 3
The start and end configuration of the experiment cases.
Case
Point
Angle
Start
End
Start
End
1
(0.3, −0.8, 0.2)
(0.8, 0.1, 0.3)
−π/2
0
2
(0.8, 0.1, 0.3)
(0.1, 0.9, 0.2)
0
π/2
3
(0.1, 0.9, 0.2)
(0.5, 0.6, 0.9)
π/2
0Figure 17
Motion path of the experiment cases.The generated sloshing, as measured by level sensor is shown in Figure18. It is measured for ten seconds from the start of motion. The motion time, residual vibration, and vibration reduction of the three cases are shown in Table 4. The motion time is the total time required to move from the start point to end point. The maximum residual sloshing for both unshaped and ZVD shaped trajectories is measured from the ZVD shaped motion time until the end of measurement, so as to make fair comparison.Table 4
Result of experiments.
Case
Motion time[s]
Max residual sloshing[mm]
Sloshing reduction[%]
Unshaped
ZVD
Unshaped
ZVD
1
1.791
2.197
25.735
6.123
76.207
2
1.783
2.189
22.1
5.929
73.172
3
1.156
1.562
45.997
9.921
78.431Figure 18
Sloshing measurement of (a) Case  1, (b) Case  2, and (c) Case  3.
(a)
Case  1
(b)
Case  2
(c)
Case  3In general, we can see that the shaped trajectories generate much less sloshing at expense of longer motion times. The additional motion time is 0.406 seconds for all cases. That amount is equivalent to the length of the ZVD input shaper, which equals one vibration period. The sloshing reduction is more or less the same in all cases, where the largest is as much as 78.431% in Case  3. The motion path of Case  3 is mostly in upward direction. This shows that the solution does not have problem in handling trajectory that contains vertical motion.


## 6. Conclusion and Further Work

This paper has examined the case of liquid container transfer system using multi joint robot arm where quick motion and minimal sloshing are concurrently addressed. The proposed framework integrates trajectory planning in three-dimensional space and input shaping for sloshing suppression. The trajectory planning algorithm, based on cubic spline optimization, has been built to generate trajectories for translation motion in three-dimensional space with joint kinematic constraints and obstacle avoidance. The sloshing suppression is built into the algorithm, based on principle of input shaping and decoupling of the robot motion into translation and hand rotation subsystems. It has been shown using simulations and experiments that the generated trajectories induce much less sloshing, thus demonstrating the benefit of the suppression strategy. Further work will mainly be directed to improve the current framework by incorporating the dynamics consideration into the algorithm. Another planned improvement is to provide redundancy by utilizing one additional joint of the robot arm for better maneuverability.



---

*Source: 279460-2015-05-18.xml*