Real-time stiffness compensation and force control of cooperating robots in robot-based double sided incremental sheet forming

In robot-based incremental sheet forming, the forming robot is displaced due to the forming forces and the comparable low stiffness. As the forming forces can not be predicted precisely, the stiffness needs to be compensated based on the measurement of a force torque sensor. While previous approaches used precalculated lookup tables, this publication presents a multi body system robot model that can calculate the displacement of the tool center point in real-time. In incremental sheet forming, the supporting robot is typically force controlled to bring superimposed stress into the forming zone. Unfortunately, this could lead to oscillations in the stiffness compensation. The presented force control approach takes the stiffness compensation into account to ensure a smooth movement of the forming robot. In a series of 30 forming experiments the effectiveness of the developed stiffness compensation and force control is validated.


Introduction
Sheet metal parts are omnipresent in mechanical engineering. While high batch sizes can be effectively produced with processes like deep drawing, the production of low batch sizes and prototypes suffers from high costs and delivery times for tools. One process that enables the economic production of sheet metal parts in small lot sizes is incremental sheet forming (ISF). In ISF, the part is formed by the movement of a universally applicable, typically hemispherical tool [1]. The tool moves along a usually spirally toolpath into a clamped sheet to form the part. This basic process variant is called single point incremental forming (SPIF) [2]. As already formed parts of the geometry are subsequently deformed by torque induced through forming forces, the geometric accuracy of SPIF remains low. One possibility to improve the accuracy of the formed part is the process variant named two point incremental forming (TPIF) [3]. By forming the part into a part dependent die, the subsequent deformation is greatly reduced. While this leads to comparable better forming results, the need for additional tooling contradicts the initial motivation of ISF. Therefore, double sided incremental forming (DSIF) has been developed which exhibits the advantages of TPIF without the need for a die by integrating a second tool into the forming process [4,5]. This supporting tool mirrors the forming tool on the opposite side of the sheet, locally building a die (see Fig. 1). This way, forming results comparable to TPIF can be achieved. Furthermore, the use of a second, identical tool enables a switch of the tools' tasks during the forming process. The possibility of forming in two different directions can be utilized to form parts that are not formable with SPIF making it a highly flexible ISF variant [7].

DSIF parameters
The use of a second tool in DSIF is accompanied by two unique process parameters, not present in SPIF and TPIF. During the forming process, the position of the forming tool is dependent on the part's geometry. However, the position of the supporting tool can be adjusted (see Fig. 2). It is required that the connection line of both tool center points (TCP) is identical to the sheet normal vector to not alter the geometry. Therefore, the supporting tool can be rotated around the forming tool, specified by the support angle . The support angle starts in the position 0 , where both tools stand directly opposed to each other and the longitudinal tool axes are in alignment. By rotating the supporting tool around the forming tool, the support angle increases. The last point where the sheet normal stops passing through the tool center points is specified as max . The position of the point depends on the wall angle which lies between the clamping frame plane and the side surface of the part. Therefore, the support angle is given in dependence of the wall angle ranging from Fig. 1 Process principle of DSIF. While the forming tool actively forms the sheet, the supporting tool locally builds a die reducing the subsequent deformation. Adapted from [6] Fig. 2 Function scheme of the support angle . While the position of the forming tool is fixed, the supporting tool can be rotated until max is reached. As this position depends on the wall angle , the supporting angle is specified in percent of max respectively . Adapted from [6] 1 3 0 = 0 % to max = = 100 % . During the forming process, the forming tool is position controlled, shaping the part. While the supporting tool can also be position controlled, it is typically force controlled, yielding various advantages. By applying the support force F S towards the forming tool, the contact to the sheet is ensured. Otherwise, due to position inaccuracies of the forming setup, loss of contact of the supporting tool can occur resulting in springback of the sheet and marks on the surface. Furthermore, the combination of support force and support angle induces superimposed stress into the forming zone [8,9]. This leads to a purposeful increase of forming accuracy and formable wall angle [10,11].

Robot based ISF
In ISF, the forming tool is typically mounted on a NC machine benefiting of its high positioning accuracy and stiffness [3]. While these are important advantages, the integration and synchronization of a supporting tool presupposes complex adjustments of the setup. Therefore, the application of NC machines is normally limited to SPIF although there are custom built DSIF NC machine setups [5,12]. An alternative is robot based ISF (Roboforming) where the tools are mounted on industrial robots. The use of industrial robots brings the advantages of a large working area and an already wide distribution. This reduces the required investment costs for ISF, making it attractive for industrial use [4]. Furthermore, industrial robots have built-in methods for tool synchronization and sensor integration. While NC machines exhibit a high structural stiffness, industrial robots have a comparably low structural stiffness because of the open kinematic chain [13]. This leads to positioning inaccuracies of the forming tool during the forming process, resulting in a deviation of the formed part. Therefore, industrial robots of the heavy load class are used as they have the highest structural stiffness. Nevertheless, even heavy load industrial robots show positioning inaccuracies due to the forming forces. Laurischkat measured a displacement of the forming tool of over 3 mm at a typical forming force of around 2000 N using two KUKA KR360 robots [14]. To ensure a precise forming, the stiffness of the industrial robots needs to be compensated. The Roboforming setup at the Chair of Production Systems (Lehrstuhl für Produktionssysteme, LPS), that is investigated in this publication, consists of two KUKA KR600 industrial robots and two KUKA KRC4 controllers (see Fig. 3). The use of a modern robot controller enables the integration of more advanced stiffness compensation algorithms to further increase the geometric accuracy of the process. This publication presents an online stiffness compensation algorithm for both industrial robots in combination with force control of the supporting tool.

Robot stiffness compensation
Load dependent positioning inaccuracies of industrial robots mainly result from joints as links are distinctly more stiff. Due to the dependency of the positioning accuracy on the joints, the displacement of the tool differs for various poses [15]. Therefore, the stiffness compensation needs to address the compliances of the individual joints. In literature, this is typically divided in modeling of static and dynamic stiffness as stiffness compensation is usually applied to robots used in milling. While milling forces are around 10 times lower than forming forces in ISF, they exhibit high frequencies, making consideration of dynamics necessary [16]. However, these can be neglected in ISF as the forming forces show no rapid changes when using low tool speeds [17]. Therefore, stiffness compensation in ISF focuses on static displacements due to the forming and supporting forces. The possibilities to simulate the ISF processes and the resulting forces are still limited. Li et al. observed an average error of 11% in terms of predicting the forming forces of a SPIF process for a conical part [18]. For more complex parts and the more difficult to simulate DSIF process these prediction errors can be expected to much bigger [19,20]. Therefore, the robot stiffness needs to be compensated online based on the measurement of a force torque sensor (FTS). While many approaches for stiffness compensation focus on specific milling poses [21], in ISF a universal approach is required to keep the process' flexibility. Abele et al. implemented a flexible joint multi body system robot model for stiffness compensation in ISF [22]. The model calculates the acting torques in each joint. Based on beforehand measured curves, describing the angular displacement of each joint depending on the torque, the displacement of the TCP is calculated. As the model was too computationally heavy to execute in real-time, lookup tables for a wide amount of force and pose combinations were calculated. Afterwards, an online stiffness compensation based on the lookup tables and the measured forming forces was implemented and validated [14]. This publication presents an improved version of the flexible joint multi body system robot model which is capable of real-time execution to further improve the robot accuracy.

Combined stiffness compensation and force control
Concurrently to the stiffness compensation of both robots in ISF, the supporting robot is force controlled (see Sect. 1.1). During the process, stiffness compensation and supporting force act in opposite directions. As both have been implemented as individual control loops before, they directly influence each others behavior, leading to problems. Thyssen observed that a forming experiment with a desired supporting force of 200 N led to standard deviation of 43.19 N although a mean of 200.02 N could be achieved [23]. Therefore, in this new approach both control loops are connected to prevent oscillation. The developed external real-time controller is implemented using Simulink Real-Time, executed on a standard industrial computer. A schematic of the signal flow is displayed in Fig. 4. For data transfer with the robot controllers, KUKA robot sensor interface (RSI) with an interpolation cycle of 4 ms is used. RSI enables the real-time adjustment of cartesian position and axis values based on sensor data or external signals. The adjustments are actuated during the following interpolation cycle. Via user defined protocol (UDP), an Ethernet connection between the external controller and both robot controllers is established. Both robot controllers send their internal cartesian position, all six axis positions and the force and torque values measured by the FTS to the external controller. Based on these values, a multi body system robot model calculates the displacement and actual position of both robots. In case of the forming robot, a position controller calculates the position offset that is sent to the robot controller and actuated in the next interpolation cycle. To prevent oscillation, the position and force controller of the supporting robot also takes the control values of the forming robot into account. Furthermore, the desired supporting force and its direction, which were calculated beforehand during the process and toolpath planning, are used to calculate the position offset. Afterwards, the position offset is sent to the robot controller. During the forming process itself, the robots follow a spiral toolpath calculated for the finishing process step of a 5-axis milling machine. The individual toolpath points are connected with linear movements. Their distance is altered to keep the approximation error in radii under 0.1 mm . As the industrial robots smooth their movement, the path error is even lower. The movement of both robots is synchronized by using the software package KUKA RoboTeam. This ensures that both robots reach each individual toolpath point at the same time.

Joint stiffness modeling and measurement
The later established stiffness compensation model is, in principle, a rigid link flexible joint model [24]. Although the links are considered rigid, their stiffness is also indirectly taken into account. This is done by assuming a total stiffness c * for all joints numbered by i. This total stiffness c * is composed of the joints individual stiffness and the structural stiffness, denoted by s, of the following link [14]: This offers the advantage of a rather low model complexity while also considering both, joints and links. Furthermore, c * can be regarded in different spatial directions (see Fig. 5). The joint's intended rotational direction r is determined by the joint's rotational axis. Orthogonal to the rotational axis, two tilting axes with the tilting directions t1 and t2, also denoted tj, are assumed. The appliance of a torque will lead to a displacement of the joint in the torque's direction according to the corresponding stiffness. Therefore, to calculate the torque and direction dependent stiffnesses of the joints, the torque and direction dependent angular displacement needed to be measured during experiments. This was done by applying a defined force at the flange, leading to torques at the individual joints (see Fig. 6). For applying the force, a pneumatic cylinder was mounted at the robot's flange in different orientations. The measured robot had its drives in control, staying in position. However, it is worth noting that no difference to test experiments with activated brakes could be observed. The pneumatic cylinder exhibits a piston with a push pull disc that stayed inside a piston intake. The piston intake was mounted on the flange of the second robot that had its brakes activated. By increasing or decreasing the pressure inside the cylinder, defined forces in two directions could be applied as the piston pushed or pulled at the intake. During the measurement, the measured robot was posed so that leverage maximized for the measured joint. This enabled the appliance of high torques. For the given example in rotational direction of joint 2 (see Fig. 6) the leverage would be the length of link 3 l L3 . Furthermore, the following joints 3-6 were posed so that their leverage was zero to make their influence negligible. However, this was not possible for previous joint 1 so that it would also be displaced which was taken into account during the measurement. The rotational displacement ,r2 of joint 2 was measured by the optical coordinate measurement system Metris K600. A measurement LED was glued on end of link 3 which was tracked during the appliance of the force. For the measurement of shorter links as link 6, extensions were used to reduce the influence of measurement noise of the LED tracking. Furthermore, a coordinate system (CS) was defined in the center of joint 2 using a space probe. This coordinate system was coupled with a couple of LED glued on link 2. By measuring the displacement of the measurement LED in the coordinate system represented by the CS LED, the influence of previous joint, in this case joint 1, is eliminated. This enables the precise measurement of ,r2 . As the robots are equipped with FTS, the applied force was measured and could be converted to torque. Therefore, a curve representing the correlation of angular displacement to torque ( ) could be created for every joint and direction (see Fig. 7). During the one hour long measurement, the pressure has been increased slowly and steadily, then decreased and afterwards increased again Signal flow between the two robot controllers and the external controller executed on an industrial computer to implement the stiffness compensation and force control. The external controller sends position adjustments via the Ethernet UDP connection to both robot controllers. These are processed by the RSI until the force reached the maximum payload in both directions, and finally finished at zero. Meanwhile, the applied force was measured by the FTS and the LEDs were tracked by the optical coordinate measurement system. Afterwards, the data was synchronized. Due to the increase and decrease of the force, this yielded the hysteresis curve of the angular displacement. However, to keep the control expenditure low, polynomial regression was used to create the curves in Fig. 7, which are used in the developed stiffness compensation approach. Although, the collected data will enable the modeling of hysteresis effects of the joints in the future. By investigating the curves, it can be seen that no rotational stiffness was measured for joint 6. As the tools used in the forming process are rotationally beared, this has no influence on the process. Furthermore, the curves for both tilting directions are identical. While this is comprehensible regarding the joint's influence on the stiffness, as it is symmetrical, the structural stiffness should be different at first sight. However, this underlines that the structural stiffness of the heavy load robots can be considered rigid.

Multi body system robot model
After determining the stiffness of all joints of the robots, a physical multi body system model has been build up to predict the displacement of the TCP due to the forming forces. Matlab Simscape Multibody has been used for this purpose as it can be integrated directly into Simulink. In combination with Simulink Desktop Real-Time, this enables calculating the displacement in real-time while also controlling the system. In Simscape, the robot was created by links and joints according to the specifications of the manufacturer (see Fig. 8). Links were created rigid with the base being grounded. A world coordinate system was created in the center of it. By utilizing the base and tool Depiction of the torque dependent total stiffness c * of the joint i in rotational direction r and the two orthogonal tilting directions t1 and t2. Adapted from [14] Fig. 6 Exemplary measurement of the angular displacement ,r2 of joint 2 as a result of the torque caused by the applied force and the leverage l L3 . By measuring the LED and the applied force the stiffness of joint 2 can be determined. Adapted from [14] data of the robot controller, all coordinates can be transformed into the robot coordinate system. The joints were each implemented in two different ways with a dummy link of vanishing length between them. The rigid joint is actuated with the axis values sent by the robot controller. Neglecting the second version of the joints, the multi body system model would simply calculate the same TCP position as the robot controller. However, the second joint is flexible although it can also be actuated. To calculate the displacement of the TCP, the forces measured by the FTS are added to the model. After solving the multi body system, the torques acting in each joint and direction are available. These can be to used to actively actuate the flexible joints according to the measured curves in Fig. 7. However, this will lead to problems as it is an iterative process. With every small change of the TCP position, the leverage of every joint is also altered leading to a new displacement of the TCP. As the robot controller expects a solution in 4 ms , this might be too computationally heavy. Therefore, two different versions of the Simscape model are used. The rigid version solely calculates the torque acting in each joint according to the measured forces and the pose. This will lead to small errors as influence of the displacement on the leverages is ignored. However, this can be neglected due to magnitude of the leverage, respectively the link length compared to to the displacement of the TCP. Afterwards, the joints of a second model are actuated with the angular displacement in all three directions additionally to the axis values sent by the robot controller. Therefore, the real TCP due to displacements by the forming forces is calculated. This can be used for stiffness compensation control.

Stiffness compensation control
The two multi body system robot models are implemented in a control algorithm in Simulink Desktop Real-Time to compensate the stiffness of the forming and the supporting robot during the process. The corresponding block diagram is shown in Fig. 9. The actual force and axis values, denoted by m, are measured and sent by the robot controller, building the input. The robot controller sends and expects data every 4 ms, in the interpolation cycle. However, the control algorithm is executed every 0.5 ms to guarantee an answer in time. This is especially important as peaks in the network traffic of up to 1-2 ms could be observed. Afterwards, the rigid multi body system robot model is used to calculate the angular displacement and the cartesian TCP position X m , which is similar to the value of the robot controller. However, X m could also be sent by the robot controller. As the robot controller is running at its computational limit, this is used to reduce its calculation load. Afterwards, m and are used in the flexible multi body system robot model to calculate the real and displaced robot position X. Therefore, the cartesian displacement X can be calculated: As X m is the position of the robot controller, formerly actuated position offsets X a by the RSI are already part of this value. Therefore, they need to be subtracted from the displacement to receive the position error X e . As the the compensation needs to point in the opposite direction of the displacement, this leads to: Rapid direction changes of the compensation values due to measurement noise of the FTS would lead to a corresponding movement of the robot. As the acceleration and deceleration would result in measured forces at the FTS, this could lead to an oscillation. To prevent this and smooth the movement, X e is filtered. High-pass filtering in real-time delays the signal, whereby the control algorithm became unstable. Therefore, a moving average filter (MAF) with a (4) X e = X a − X range of 320 time steps is used. While this might seem long regarding a filter range of 160 ms respectively 40 interpolation cycles of the robot, there are no rapid changes of the forming force in ISF. Afterwards, the filtered position error X e,f is processed by a PID controller, treating each cartesian direction individually. As the robots dynamics are considered as a black box, its values were determined in an iterative approach. For a PID controller equation of  Table 1 were determined. The controller's output X c could be sent to the robot controller afterwards. As the RSI would actuate the offset even when it surpasses the capabilities of the drives, this could lead to termination of the program if the offset is too big. Therefore, a self-incrementing approach has been chosen where the output of the PID controller resembles only the change of the offset in one time step instead of the total value. This way, the maximum change can be capped, which prevents overly high accelerations. Although this was only necessary during the testing phase, the PID controller is parameterized to function with the self incrementation. This offset value X a is later sent to the robot controller and actuated by the RSI.

Design of experiments
For validating the stiffness compensation, a dedicated geometry has been designed (see Fig. 10). The geometry exhibits flat, convex and concave surfaces to cover a wide range of movement conditions. Furthermore, each corner is rounded with a radius of 12 mm. This leads to a comparatively fast change in the direction of the forming forces. This way, the direction of the stiffness compensation needs to rapidly change during the process and the (5) X c = K P 1 + K I T s 1 z(X e,f ) − 1 + K D N 1 + NT s 1 z(X e,f )−1 dynamic behavior can be tested. The geometry itself is extruded with a wall angle of 50 • . To validate the stiffness compensation, an experimental series consisting of 24 individual experiments has been designed. The full design of experiments is given in Table 2. In every experiment, the same geometry has been used to calculate the robots' tool path. To validate the influence of both, the forming process itself and the stiffness compensation, four sets with six experiments each were built. In them there was either a metal sheet clamped and therefore formed or the  experiment was carried out without a metal sheet. This results in missing forming forces as the tools are only moving through the air without contact. Therefore, occurring forces were only resulting from inertia. Furthermore, the experiments were carried out with activated or deactivated stiffness compensation to compare the displacements. Additionally, step depth and supporting angle were varied. In the experiments with active forming, a DC04 metal sheet with a thickness of 0.8 mm has been formed. The tool speed was 20 mm/s in all experiments. During the experiments, the external control computer logged the target TCP positions of both robots without any offsets actuated by the RSI. These positions were determined during the earlier path planning process. To calculate the positioning error, the robots' TCP positions were tracked with the optical coordinate measurement system Metris K600 during the process. A plate with five measurement LED has been mounted on each robot flange. These LEDs were coupled with the TCP utilizing a space probe. Therefore, the measured LED positions during the process could be used to calculate the TCP positions of both robots.

Experiment evaluation
To calculate the positioning error of the tools during the experiments, the target position and the measured position were synchronized according to the measurement frequencies and the start of the movement. Unfortunately, random desynchronization of the data could be observed. Although the error remains small, up to 0.1 s, this equals 2 mm at the given movement speed. Either the controller of the coordinate measurement system or, more probably, the external controller did not keep a constant measurement frequency. As over three million data points were gathered in total, manual synchronization was not possible. Therefore, another approach was used. The measured TCP positions were used to calculate their minimum distance to the CAD model. As this shows always along a surface normal vector, displacements inside the plane of the surface can not be detected.
Although it is worth noting, these would not influence the forming result as the formed geometry would be the same if the displacement occurs in the plane of the surface. Unfortunately, this evaluation can only be done for the forming robot. Due to the incremental tool path and the supporting angle, the position of the supporting robot depends on the position of the forming robot, and not directly on the CAD model. Therefore, only the forming robot has been investigated. As the supporting robot is exposed to smaller forces, it can be taken as granted that its displacement would be smaller than that of the forming robot. The measured positioning errors of the forming robot are given in Table 3. For E01-E06, where no sheet was formed and the stiffness compensation was deactivated, a mean error of close to 0 mm could be observed. The mean average error (MAE) and the standard deviation are, being around 0.2 mm, typical for industrial robots of this size. However, the minimum and maximum are rather large, around ±0.5 mm . Although the path accuracy is much worse, this seems somewhat high for the low tool speed. Especially, as the robots are absolutely calibrated. The reason for these errors is the heavy computational load of the process. As the robot controllers need to synchronize their movement, send data to the external controller and alter the movement according to the RSI (although this applies only for E07-E12 and E19-E24), the robot controllers operate at their computational limit. This is amplified by the used tool paths. The robot tool paths consist of linear movements, which results in a small distance between the tool path points at smaller radii. As a result of the computational load, the KUKA smartPADs became unusable during the process and for almost half an hour afterwards as the screens stopped refreshing and button presses were not recognized. Furthermore, it occurred in all 30-90 min long experiments around five times that the movement stopped for a short amount of time as the robot controllers could not finish their path planning in time. This deceleration and acceleration caused the minimum and maximum errors. As the stiffness compensation control algorithm features self incrementing elements, a stop of the movement would lead to wrong calculations. Therefore, the self incrementing elements were coupled to the movement and only executed if the robot was moving. By comparing the positioning errors of E01-E06 to E07-E12, where the stiffness compensation was activated although no sheet was formed, a slight improvement can be observed. This seems surprising as no stiffness compensation is needed due to the nonexistent forming forces. The influence of the inertia can be neglected because of the low tool speeds. However, the plate, on which the measurement LED are glued, was not registered in the load data. Therefore, the compensation of its load is causing the improvement. What also can be observed, is a worsening of the positioning accuracy in E12. For unknown reasons a higher amount of stops occurred due to computational overload in this experiment. As this is no result of the stiffness compensation algorithm itself, its applicability for none and low forces can be taken as granted.
In E13-E18 metal sheets were formed without compensating the robot stiffness. Due to the forming forces of over 1000 N, huge displacements and therefore positioning errors were expected. They could be observed as the mean MAE of the six experiments is 0.4277 mm while the mean maximum error reached 1.2076 mm . Furthermore, only less than half of the path reached an absolute accuracy of less than 0.3 mm.
These huge displacement should be eliminated by activating the stiffness compensation in E19-E24. As it can be seen in Table 3, the stiffness compensation significantly reduced the positioning error, even surpassing E01-E06 due to the plate with the measurement LED. The mean values of the set experiments indicate even better results than in E07-E12. However, this is not completely the case as E12 worsens the mean value because of the high amount of stops. Looking at the individual values, it can be seen that E19-E24 overall exhibit slightly worse values than E07-E12. However, as the positioning error improved, even compared to E01-E06, the stiffness compensation is working as intended.

Force control
In ISF, the supporting robot is typically force controlled.
In combination with the supporting angle, superimposed stress can be induced, greatly improving the forming result. Force control algorithms are generally divided in explicit and implicit force control [25]. In explicit force control, the robot's position is directly altered based on the force error. This is not possible in ISF as the sheet thinning can not be predicted precise enough to calculate a new position in which the desired supporting force is reached. Furthermore, the positioning accuracy of the forming robot also has an influence, as both tools press against the metal sheet on opposing sides. Therefore, an implicit force control algorithm is used where a position error is calculated and controlled. Force control in ISF has two challenges. The direction of the force vector is not constant as the supporting force should act towards the forming tool while the supporting tool follows it on the outside of the sheet. Since the supporting force causes forces at the forming tool which are compensated by the stiffness compensation, this can lead to oscillations. These challenges are addressed by the presented force control approach. For achieving this, the distance between the tools, calculated considering their displacements, is incrementally altered based on the force error. The corresponding block diagram is given in Fig. 11. Inputs of the control system are the force F m , measured by the FTS, and the Euler angles E, which have been precalculated during the tool path planning. These rotate the coordinate system so that the Z axis points in direction of the surface normal and the Y axis towards the next tool path point. As a first step, F m is transformed into the coordinate system specified by E by utilizing a rotation matrix. The Z component of the output vector F m,t resembles the current supporting force. By subtracting the desired supporting force F S , the force error F e is calculated. The force error is afterwards multiplied with the gain factor F g to calculate the desired change of the distance between the tools in the current time step. This is added to the desired distance between the tools D F . The value of F g is given together with all other control values in Table 4.
To calculate the distance error D e , the current distance D between the two tools is needed. Therefore, the current position of the forming robot X FR and the supporting robot X SR , taking the displacement into account, and their radii R FR∕SR are used to calculate D. As D does not exhibit rapid changes, its value is filtered by a moving average filter (MAF), considering 200 time steps, to smooth the movement. Although, this step is not necessary. The filtered current distance D f is then used together with the already actuated distance D a to calculate the distance error: Afterwards, a PD controller of the same structure as already used in stiffness control (see Eq. 5) calculates the controlled distance D c . This resembles the change of the actuated distance D a . In the last step, D a is converted into a three-dimensional position offset X a that is added to the position offset of the stiffness compensation. This value is actuated by the RSI. To calculate X a , the initial coordinate transformation is inverted. For this, the vector (0, 0, D a ) and the Euler angles are used. As the Euler angles were precalculated for every tool path point, they only change when a new tool path point Fig. 11 Block diagram of the developed force control algorithm. To achieve a constant supporting force, the distance D between the tools is altered based on the force error. The distance offset D a is trans-formed into a three-dimensional offset X a that is sent to the robot controller, added to the offset of the stiffness compensation and actuated by the RSI is reached. This rapid change of the coordinate system would lead to jumps of X a and a high acceleration of the robot. Therefore, E is filtered by a moving average filter with a time span of 200 time steps.

Design of experiments
For validating the force control algorithm, E19-E24 were repeated as the supporting force requires an actively formed sheet metal. Otherwise, the two tools would directly press against each other what might cause damage to tools and robots. Furthermore, the stiffness compensation was activated for both robots to resemble the normal process conditions. In the validation experiments E25-E30 (see Table 5) a constant supporting force of 200 N was used as this usually yields the best forming results [26]. During the forming process, the supporting force increases the forming force as both forces act in different directions. Therefore, the displacement of the forming tool will also be increased, especially if it oscillates. To measure the influence of the supporting force on the positioning accuracy of the forming robot, the position of the forming tool was tracked by the optical coordinate measurement system (Table 6).

Experiment evaluation
During the forming experiments E25-E30, the actuated supporting force has been measured. When looking at the statistical evaluation, it can be seen that in all experiments a mean supporting force of close to the desired 200 N has been achieved (see Table 7 [23]. This was achieved using the software package KUKA Force-TorqueControl by forming a simpler geometry with less rapid direction changes. Comparing this to the validation experiments, it can be seen that the control quality of the supporting force mainly depends on the supporting angle. While similar results were achieved for a supporting angle of 0 % , the standard deviation improved for a supporting angle of 50 % . However, at a supporting angle of 100 % the standard deviation was significantly larger. The step depth has a negligible influence on the control quality of supporting force as there are only small differences between the experiments with the same supporting angle. One exception is the maximum supporting force in E25 that reached over 1000 N. This results from the aforementioned stopping behavior due to computational overload. Nevertheless, E29 exhibits a supporting force of over 800 N which is clearly too high. By looking at the minimum supporting force value, it can be seen that in all experiments the supporting tool was close to lose or even lost the contact to the sheet metal.
The negative values for the supporting force are physically impossible and result from white noise and filtering that delayed the coordinate transformation of the force measured by the FTS. The reasons for these large supporting force errors can be derived from Fig. 12. In every experiment the measured support force exhibits really good results on the two flat surfaces and on the convex face. However, the supporting force error is large in the small radii and on the concave surface. This applies especially for E29 and E30 which matches the results in Table 7. In every experiment a small increase of the supporting force can be observed at every radius. This is expected as the coordinate systems of both robots will not be perfectly aligned. Although they were measured using an optical coordinate measurement system, an offset of approximately 0.1-0.2 mm can not be prevented as the positioning accuracy of the robots has an influence on the position of the coordinate system. When the robots rapidly change their moving direction in the radii, the offset of coordinate systems will lead to an increase or a decrease of the supporting force. The reaction of the controller could lead to oscillation which applies especially for the concave surface in Fig. 12a, b. In the other experiments, much less oscillation can be noticed. Although, Fig. 12e, f exhibit huge peaks. This could be explained by a bad alignment of the coordinate systems in the XY plane as this would have the greatest impact at supporting angle of 100 % when the tools are next to each other, instead of opposed to each other. Nevertheless, the bad control behavior can not be solely explained by an imperfect alignment of the coordinate systems as the controller should be able to minimize the error. As the developed force control algorithm does not consider the dynamics of the robot, it does not respond fast enough on the error or oscillates itself as in E25 and E26 in the concave radius. This is amplified by the calculation of the supporting force direction. The direction is precalculated for every tool path point and filtered to prevent jumps as described in Sect. 6. This brings a small delay into the system resulting in the supporting force being measured in a tilted direction as the position offset of the force controller is actuated. Future research will incorporate the dynamics of the robots into the force controller and improve the calculation of the support vector direction.
As the supporting force acts in direction of the forming tool, a portion of it will enlarge the forming force and therefore the displacement that needs to be compensated. Especially the occasional peaks of the supporting force that could be observed during the validation experiments will have an influence on the positioning accuracy. To determine the positioning accuracy, the forming robot has been tracked by the optical coordinate measurement system during the experiments and the gathered position data has been evaluated (see Fig. 6). Compared to the experiments E19-E24 without force control, no remarkable differences of the positioning accuracy can be observed. Although, E30 worsened a bit compared to E24 due to the supporting force peaks, but is still comparable to E06 where no sheet was formed and no compensation was active. Furthermore, the mean error slightly worsened by around 0.05 mm in all validation experiments. While this might result from the supporting force, it is more probable that a slightly different calibration of the coordinate system of the optical coordinate measurement system is the reason, as the MAE is nearly unchanged. Therefore, the developed stiffness compensation is able to ensure a high positioning accuracy of the forming tool, independent of the acting supporting force.

Conclusion
The forming forces in ISF lead to a displacement of the used heavy load industrial robots due to their comparably low stiffness that could be demonstrated in an experimental series. To compensate for the stiffness, a multi body system robot model has been developed which can calculate the angular displacement of the individual joints and therefore the displacement of the TCP. The model is capable of realtime execution in under 4 ms on an external controller. In a series of experiments, the effectiveness of the approach could be validated as the positioning accuracy of the forming robot during the forming process reached and even surpassed the positioning accuracy when no sheet metal was formed and therefore no forming forces were acting. As superimposed stress increases the geometric accuracy of the formed parts, the supporting robot is force controlled to apply a force on the backside of the sheet. This supporting force acts against the stiffness compensation of the forming robot which could lead to oscillations. This is addressed in the developed force control approach. The quality of the force control algorithm has been validated in experiments and mostly even surpasses the built-in solution of the robot only every 300th data point is displayed. The differences in the outer dimensions result from the different tool paths due to the different supporting angles. The tool path runs counter clockwise around the geometry manufacturer. However, rapid changes of the movement direction lead to peaks of the supporting force. These result from not considering the robot dynamics and delayed calculation of the supporting force direction and will be addressed in future research. Nevertheless, the positioning accuracy of the forming robot stays unchanged despite the supporting force peaks clarifying the high quality of the developed stiffness compensation approach.