## 4 Numerical Simulation  
### System Parameters  
In this simulation analysis, we only consider the flexibility of the first three joints of the robot because the moment arms of the cutting force about the last three joints are much shorter than those about the first three joints. The aforementioned coordinate systems in this simulation are described as follows. The inertial coordinate system \(\{o_xyz\}\) and the local coordinate systems \(\{o_{i}\_x_{i}y_{i}z_{i}\}\) (\(i=1,\,2,\,3\)) of the robot are established as shown in Figure 5a. The inertial coordinate system is located at the base of the robot with the positive \(z\)-axis vertically upward and the \(y\)-axis perpendicular to the plane containing the second and third links. The \(z_{i}\)-axis of the local coordinate systems are along the rotation directions of the corresponding joints; the axis \(x_{1}\) is parallel to the axis \(x\) when the robot is stationary, and the \(x_{2}\) and \(x_{3}\) axes are along the link length directions. The directions of the other reference axes of \(\{o_{i}\_x_{i}y_{i}z_{i}\}\) (\(i=1,\,2,\,3\)) are determined from the right-hand rule. The local coordinate system of the work table \(\{o_{c}\_x_{c}y_{c}z_{c}\}\) is established at the middle position near the robot, with the \(x_{c}\) and \(z_{c}\) axes parallel to the \(x\) and \(z\) axes respectively and the \(y_{c}\) axis is determined from the right-hand rule. The workpiece and the milling path are shown in Figure 5b. The workpiece coordinate system \(\{o_{w}\_x_{w}y_{w}z_{w}\}\) is established at the bottom of the workpiece, with the \(x_{w}\) and \(y_{w}\) axes parallel to the length and width and also parallel to the \(x_{c}\) and \(y_{c}\) axes respectively and the \(z_{w}\) axis is determined from the right-hand rule. Figure 4: The calculation flow chart of \(d_{m}\). The geometric, inertial, and structural parameters of the robot are described as follows. The link length vectors of \(o_{10}\)\({}_{2}\), \(o_{2}\)\(o_{3}\), and \(o_{3}\)\(o_{t}\) can be written in local coordinate systems of robot as \(\mathbf{I}_{12}=[0,0,l_{1}]^{T}\), \(\mathbf{I}_{23}=[l_{2},\ 0,\ 0]^{T}\), and \(\mathbf{I}_{12}=[l_{3},\ 0,\ 0]^{T}\), where the link lengths are \(l_{1}\) = 0.2 m and \(l_{2}=l_{3}\) = 0.8 m. The masses, position vectors of the centers of gravity (CG), and inertia tensors of the links are listed in Table 1. Note that these vectors and tensors are written in the corresponding local coordinate systems. The torsional stiffnesses and damping of the three joints can be obtained by a stiffness identification experiment [45; 46]. However, our system is still being built, and thus, the joint stiffness values available for a robot with similar payload and size are adopted, i.e., \(k_{1}\) = 2.4 \(\times\) 10\({}^{5}\) Nm/rad and \(k_{2}\) = \(k_{3}\) = 10\({}^{5}\) Nm/rad. A damping ratio \(\zeta\) of 0.06 is used to compute the damping matrix based on the modal test data for industrial robots given in Reference [46]. Then, the ranges of the joint angles are \(\theta_{1}\in[-\pi,\pi]\), \(\theta_{2}\in[-\pi/2,\pi/2]\), and \(\theta_{3}\in[-\pi/3,\pi/3]\). The geometric parameters of the work table and the workpiece are described as follows. The position vector of the origin of the work table coordinate system \(o_{c}\) are written as \(\mathbf{p}_{c}\)= \([p_{x},\ p_{y},\ p_{z}]^{T}\), where \(p_{x}\), \(p_{y}\), and \(p_{z}\) are equal to 0.6 m, 0 m, and 0.5 m respectively. The length, width, and height of the work table are 1.16 m, 0.8 m, and 0.5 m, respectively, and the inclination angle of the fixture \(\theta_{p}\) can change from 0 to \(\pi/3\). The length \(L\) and width \(D\) of the milling path are 0.16 m and 0.08 m, respectively. We select 20 intermediate points uniformly distributed along the milling groove.