# Optimal Energy Consumption for Mobile Manipulators Executing Door-Opening Task

Authors: Changyou Ma; Haibo Gao; Liang Ding; Jianguo Tao; Kerui Xia; Haitao Yu; Zongquan Deng
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

As a substitute for humans, the mobile manipulator has become increasingly vital for on-site rescues at Nuclear Power Plants (NPPs) in recent years. The high energy efficiency of the mobile manipulator when executing specific rescue tasks is of great importance for the mobile manipulator. This paper focuses on the energy consumption of a robot executing the door-opening task, in a scenario mimicking an NPP rescue. We present an energy consumption optimization scheme to determine the optimal base position and joint motion of the manipulator. We developed a two-step procedure to solve the optimization problem, taking the quadric terms of the joint torques as the objective function. Firstly, the rotational motion of the door is parameterized by using piecewise fifth-order polynomials, and the parameters of the polynomials are optimized by minimizing the joint torques at the specified base position using the Quasi-Newton method. Second, the global optimal movement of the manipulator for executing the door-opening task is acquired by means of searching a grid for feasible base positions. Comprehensive door-opening experiments using a mobile manipulator platform were conducted. The effectiveness of the proposed method has been demonstrated by the results of physical experiments.


---

## Body

## 1. Introduction

Mobile manipulators, as a replacement for humans, play a key role when they perform rescue tasks in the extreme environment of a nuclear power plant (NPP), such as door-opening and turning a valve. The energy optimization of rescue tasks performed by a robotic arm, in particular that of the fundamental task of door-opening, is of great importance when the energy supply is restricted by the capacity of the battery with which the robot can be equipped, unlike in traditional methods where energy is supplied via a cable.The Fukushima Daiichi NPP accident in Japan on March 11, 2011, was triggered by an earthquake of magnitude 9.0 and the resultant tsunami marked the beginning of the worst nuclear accident of the last two decades [1–3]. Nuclear power is an important resource and (NPPs) have recently undergone rapid development in China [4]. The complexity of the NPP’s structure, the radioactive working environment, and some special characteristics, such as high temperatures and high pressure, render nuclear disaster rescue very difficult. The dangers caused by such an accident mean that mobile robot plays an important role in the NPP rescue process. After the Fukushima Daiichi NPP accident, the Defense Advanced Research Projects Agency (DARPA) initiated a new challenge in 2012, the DARPA Robotics Challenge (DRC) [5], including eight rescue tasks to test the capabilities of teams, as well as of individuals. The tasks were ranked by DARPA in terms of difficulty as Valve (easiest), Terrain and Hose (easier), Door, Debris, Wall, and Ladder (harder), and Vehicle (hardest) [6].In this paper, we focus on the energy consumption of a mobile manipulator performing the door-opening task. Nagatani and Yuto proposed a method that allows a mobile manipulator, named the “YAMABICO-10” robot, to open a door and pass through the doorway. In their approach, they applied to the mobile manipulator control system the concept of action primitives, which control the “YAMABICO-10” robot according to sequences of planned motion primitives. However, each action primitive was designed with an error adjustment mechanism to handle the accumulated positioning error of the mobile base [7–9]. Peterson et al. presented the design and implementation of a door-opening controller using a hybrid dynamic system model; when using this simple controller, the radius and center of rotation of the door are estimated online. The results of their experiments demonstrated that off-the-shelf algorithms for force/torque control are very effective for solving the task of grasping the handle and opening the door [10]. Chung et al. [11] proposed a control strategy for the door-opening procedure executed by a service robot, PSR1, which utilized a three-fingered robot hand for grasping the door handle. Two active-sensing strategies were proposed to estimate the kinematic parameters in a real environment. An integrated strategy of motion coordination was presented based on the components of three subsystems: a robotic hand, a robotic arm, and a mobile robot. The force and position control were successfully achieved by using the contact force of the three-fingered robotic hand during the door-opening procedures.Ahmad et al. designed a modular and reconfigurable robot (MRR) mounted on a wheeled mobile platform [12, 13]. They proposed a new method that utilizes the multiple working modes of the MRR modules to prevent the occurrence of large internal forces that arise because of positioning errors or imprecise modeling of the robot or its environments. By selectively switching the joints of the MRR to work in passive mode during the door-opening operation, the controller design was significantly simplified. Zhang et al. [14] presented a multiple mode control system of a two-degree-of-freedom (DOFs) compact wrist that can work in active mode with position or torque control, or in passive mode with wrist-environment interactive force compensation. They verified in their door-opening experiments that the wrist could move freely without generating excessive internal force. Kobayashi et al. [15] designed a rescue robot series named UMRS. The robot had a special end-effector, which can grip and rotate cylindrical type and lever type door handles. The robot equipped with a door-opening system was capable of moving freely through rooms, even if there were doors between them. Klingbeil et al. [16] proposed a method that used vision to identify a small number of key positions, such as the axis of rotation of the door handle, and the end-point of the door handle, to allow a manipulator to open various types of doors without prior knowledge of their parameters. Karayiannidis et al. [17] developed an algorithm that can be implemented in a velocity-controlled manipulator, the end- effector which is equipped with force sensing capabilities. The method consists of a velocity controller, which uses force measurements and an estimation of the radial direction based on adaptive estimates of the position of the door hinge. The control action can be decomposed into an estimated radial and tangential direction following the concept of hybrid force/motion control. Endres et al. [18] presented an approach for learning a dynamic model of a door from sensor observations and utilizing it for effectively swinging the door open to the desired angle. The learned models enable the realization of dynamic door-opening strategies and reduce the complexity of the door-opening task.As proposed in [19], a power efficiency estimation-based health monitoring and fault detection method has been developed for a modular and reconfigurable robot (MRR). The power efficiency of each of the robot’s joints is measured using sensors. Luo et al. [20] presented the Lagrange interpolation method to express each joint trajectory function to realize trajectory planning that achieves energy minimization of industrial robotic manipulators. Field and Stepanenko presented an iterative dynamic programming method that is modified to perform a series of dynamic programming, passing over a small reconfigurable grid that covers only a portion of the solution space at any one pass, to plan minimum energy consumption trajectories for robotic manipulators [21]. Liu et al. [22] proposed the fourth-order Runge-Kutta method, multiple shooting methods, and traversing method to solve optimal energy trajectory planning for palletizing robot.The objective of this study was to develop an energy consumption optimization method for the door-opening procedure. The main contributions of this paper are summarized as follows. An energy consumption optimization scheme that finds the optimal base position and joint motion of the manipulator for optimizing energy consumption is presented. Since the end-effector trajectory is assigned according to the trajectory of the door handle, the problem of searching the optimal manipulator movement by using the parameters of the door’s fifth-order splines transforms to a simple parametric optimization problem at each base position. The global optimal energy-efficient movement of the robotic arm is finally obtained by exhaustively searching the entire base position grid.The rest of the paper is organized as follows. In Section2, we introduce the door-opening task and establish the energy consumption optimization objective function. Section 3 provides a description of the mobile modular robot and introduces the kinematic and the dynamics model of the manipulators. In the Section 4, we address the optimization method and describe numerical simulations. The results of our experiments are discussed in Section 5, and conclusions are presented in Section 6.


## 2. Formulation of Optimization Problem

### 2.1. Door-Opening Operation

In this section, we propose a door- opening method. The following assumptions were made:(1) the door-opening and door handle rotating directions are known (the door is opened toward the left-side and the door handle is rotated toward the right-side); (2) the door axis of rotation is perpendicular to the floor; (3) the door moves in the horizontal xy plane; (4) the mobile platform travels on the ground, which can always be adjusted for a structured laboratory environment; and (5) the axis of rotation of the first manipulator module is perpendicular to the ground.A brief explanation of the door-opening procedure is as follows.(a)
The mobile manipulator moves such that it is positioned in front of the door and grasps the door handle;(b)
The mobile base remains static. The manipulator realizes turning the door handle by tracking a planned trajectory;(c)
The mobile base still remains static. The manipulator realizes pulling the door by tracking a planned trajectory in the horizontalxy plane.A model of the door-opening procedure using the mobile manipulator is depicted in Figure1. The figure shows the origin of the reference frame, Od, set at the intersection point of the door hinge, the horizontal plane that crosses the origin of the reference frame of the mobile manipulator, Ob, and the reference frame, Oe, located at the position of the manipulator end-effector, which is used to grasp the door handle.Figure 1
Model of the door-opening procedure with a mobile manipulator.In order to plan the path of the mobile manipulator, we first need to know the accurate value of the door handle radius(rk) and the door radius (rd), the initial base position of the mobile manipulator (xb,yb,zb), and the position of the end-effector, which holds the door handle firmly (xe,ye,ze). The door motion is conformed to follow the door trajectory in the xy plane with the center of rotation at (xd,yd) and a radius rd as shown in the Figure 1. The trajectory radius of pulling the door is derived as follows:(1)xet-xd2+yet-yd2=rd2.


### 2.2. Path Planning of Door-Opening

In this section, we focus on pulling the door handle to open the door. During this procedure, the home position of the mobile base, the door radius, and the height of the door handle are measured. These measured parameters are then used for planning the path of the mobile manipulator that allows it to open the door to the desired angle.The mobile manipulator achieves the action of pulling the door by tracking the circle arc represented by the door trajectory, as shown in Figure2. The arc of the door trajectory can be acquired by using the interpolation points P1, P2, and P3. In Figure 2, we show the optimal base position (xb,yb) during the door-opening procedure. The mobile manipulator approaches the door and stops at any position near the door. The end-effector is fixed at the same position as the door handle during the door-opening procedure, and each joint of the manipulator rotates in accordance with the motion planning. The door is pulled in the horizontal xy plane, and the initial and final angles of the door are consistent. The angles of Joints 1, 2, 3, 4, 5, and 6 are denoted as θ1, θ2, θ3, θ4, θ5, and θ6, respectively. The input torques of each joint of the manipulator are expressed as τ1, τ2, τ3, τ4, τ5, and τ6, respectively. The angles and the torques are represented in matrix forms as θ=θ1,θ2,θ3,θ4,θ5,θ6T and τ=τ1,τ2,τ3,τ4,τ5,τ6T. The angle of the door is expressed as θd, and the conditions on θd at the start time, t=0, and the end time, t=te, are written as(2)θd,θ˙dt=0=0,0,θd,θ˙dt=te=θdte,0.Figure 2
Planar model of pulling door.The positions of the end-effector at pointsP1, P2 and P3 are computed as follows:(3)θd=θdoor3n,n=1,2,3P1,2,3=pxpypz=rd1-cos⁡θdy1-rdsin⁡θdz1,where rd denotes the rotation radius of the door, θd denotes the rotation angle of the door when reaching points P1, P2, and P3, and θdoor denotes the rotation angle of the door that allows the mobile manipulator to enter the doorway.The orientations of the end-effector and the door are the same during the process of pulling the door. When the position of the end-effector remains constant with respect to the door handle, as shown in Figure1, the variation in the rotation matrix during the process denoted by Re0d(0). Re0d(0) denotes the initial rotation matrix of the end-effector after turning the door handle, which can be acquired by using forward kinematics [23]. Then, at timet, this matrix is multiplied by a rotation matrix parameterized by the rotation angle of the door θd(t):(4)Re0dt=cos⁡θdt-sin⁡θdt0sin⁡θdtcos⁡θdt0001Re0d0.


### 2.3. Objective Formulation of Energy Consumption

The energetic cost is an important metric of the energy consumption of the manipulator. We experimentally evaluated the effect of different base positions of the robot on the energy consumption of the manipulator during the door-opening procedure. The input energyEm of the motor is as follows:(5)Em=∫0teUin·Imdt,where Uin is the power supply voltage and Im is the instantaneous current of the DC motor. The total energy consumed by the robot’s actuators includes the generated mechanical power (Pmech), heat power (Pheat), and power losses (Δloss). This can be described as(6)Pin=Pmech+Pheat+Δloss.The mechanical power generated by each actuator is related to the rotation angular velocity and torque. Therefore, the total instantaneous mechanical power during the movement of the manipulator can be stated as(7)Pmech=∑j=16τj·ωj,where Pmech is the total mechanical power of the actuators, τj is the torque of each actuators in N·m, ωj is the rotation angular velocity of motors in rad/s, and j is the number of manipulator joints.The torqueτ of each actuator is related to the torque constant Km, where, i is the gear ratio of the joint, η is the efficiency of the transmission mechanism, and Im is the instantaneous current; τ is provided by [24](8)τ=Km·Im·i·η.As shown in (7), the mechanical power consists of the torque of each actuators and the rotation angular velocity of the motors. However, the angular velocity range of the joint trajectory planning of the manipulator is relatively small. Therefore, the problem of energy optimization becomes the joint torque optimization problem. For the trajectory optimization procedure, the integral of squared joint torques as the cost function is known from studies in the literature [25]. We introduce the following objective function as a standard for optimization. (9)Econς=∑j=16∫t=0teτj2dt.Here, ς denotes the optimal parameter, which is chosen as the position of mobile platform base and the trajectory of the door angle. To solve the problem of energy consumption optimization during the door-opening procedure, we determine the parameters that minimize Econ:(10)ς∗=argmin⁡Econ.


## 2.1. Door-Opening Operation

In this section, we propose a door- opening method. The following assumptions were made:(1) the door-opening and door handle rotating directions are known (the door is opened toward the left-side and the door handle is rotated toward the right-side); (2) the door axis of rotation is perpendicular to the floor; (3) the door moves in the horizontal xy plane; (4) the mobile platform travels on the ground, which can always be adjusted for a structured laboratory environment; and (5) the axis of rotation of the first manipulator module is perpendicular to the ground.A brief explanation of the door-opening procedure is as follows.(a)
The mobile manipulator moves such that it is positioned in front of the door and grasps the door handle;(b)
The mobile base remains static. The manipulator realizes turning the door handle by tracking a planned trajectory;(c)
The mobile base still remains static. The manipulator realizes pulling the door by tracking a planned trajectory in the horizontalxy plane.A model of the door-opening procedure using the mobile manipulator is depicted in Figure1. The figure shows the origin of the reference frame, Od, set at the intersection point of the door hinge, the horizontal plane that crosses the origin of the reference frame of the mobile manipulator, Ob, and the reference frame, Oe, located at the position of the manipulator end-effector, which is used to grasp the door handle.Figure 1
Model of the door-opening procedure with a mobile manipulator.In order to plan the path of the mobile manipulator, we first need to know the accurate value of the door handle radius(rk) and the door radius (rd), the initial base position of the mobile manipulator (xb,yb,zb), and the position of the end-effector, which holds the door handle firmly (xe,ye,ze). The door motion is conformed to follow the door trajectory in the xy plane with the center of rotation at (xd,yd) and a radius rd as shown in the Figure 1. The trajectory radius of pulling the door is derived as follows:(1)xet-xd2+yet-yd2=rd2.


## 2.2. Path Planning of Door-Opening

In this section, we focus on pulling the door handle to open the door. During this procedure, the home position of the mobile base, the door radius, and the height of the door handle are measured. These measured parameters are then used for planning the path of the mobile manipulator that allows it to open the door to the desired angle.The mobile manipulator achieves the action of pulling the door by tracking the circle arc represented by the door trajectory, as shown in Figure2. The arc of the door trajectory can be acquired by using the interpolation points P1, P2, and P3. In Figure 2, we show the optimal base position (xb,yb) during the door-opening procedure. The mobile manipulator approaches the door and stops at any position near the door. The end-effector is fixed at the same position as the door handle during the door-opening procedure, and each joint of the manipulator rotates in accordance with the motion planning. The door is pulled in the horizontal xy plane, and the initial and final angles of the door are consistent. The angles of Joints 1, 2, 3, 4, 5, and 6 are denoted as θ1, θ2, θ3, θ4, θ5, and θ6, respectively. The input torques of each joint of the manipulator are expressed as τ1, τ2, τ3, τ4, τ5, and τ6, respectively. The angles and the torques are represented in matrix forms as θ=θ1,θ2,θ3,θ4,θ5,θ6T and τ=τ1,τ2,τ3,τ4,τ5,τ6T. The angle of the door is expressed as θd, and the conditions on θd at the start time, t=0, and the end time, t=te, are written as(2)θd,θ˙dt=0=0,0,θd,θ˙dt=te=θdte,0.Figure 2
Planar model of pulling door.The positions of the end-effector at pointsP1, P2 and P3 are computed as follows:(3)θd=θdoor3n,n=1,2,3P1,2,3=pxpypz=rd1-cos⁡θdy1-rdsin⁡θdz1,where rd denotes the rotation radius of the door, θd denotes the rotation angle of the door when reaching points P1, P2, and P3, and θdoor denotes the rotation angle of the door that allows the mobile manipulator to enter the doorway.The orientations of the end-effector and the door are the same during the process of pulling the door. When the position of the end-effector remains constant with respect to the door handle, as shown in Figure1, the variation in the rotation matrix during the process denoted by Re0d(0). Re0d(0) denotes the initial rotation matrix of the end-effector after turning the door handle, which can be acquired by using forward kinematics [23]. Then, at timet, this matrix is multiplied by a rotation matrix parameterized by the rotation angle of the door θd(t):(4)Re0dt=cos⁡θdt-sin⁡θdt0sin⁡θdtcos⁡θdt0001Re0d0.


## 2.3. Objective Formulation of Energy Consumption

The energetic cost is an important metric of the energy consumption of the manipulator. We experimentally evaluated the effect of different base positions of the robot on the energy consumption of the manipulator during the door-opening procedure. The input energyEm of the motor is as follows:(5)Em=∫0teUin·Imdt,where Uin is the power supply voltage and Im is the instantaneous current of the DC motor. The total energy consumed by the robot’s actuators includes the generated mechanical power (Pmech), heat power (Pheat), and power losses (Δloss). This can be described as(6)Pin=Pmech+Pheat+Δloss.The mechanical power generated by each actuator is related to the rotation angular velocity and torque. Therefore, the total instantaneous mechanical power during the movement of the manipulator can be stated as(7)Pmech=∑j=16τj·ωj,where Pmech is the total mechanical power of the actuators, τj is the torque of each actuators in N·m, ωj is the rotation angular velocity of motors in rad/s, and j is the number of manipulator joints.The torqueτ of each actuator is related to the torque constant Km, where, i is the gear ratio of the joint, η is the efficiency of the transmission mechanism, and Im is the instantaneous current; τ is provided by [24](8)τ=Km·Im·i·η.As shown in (7), the mechanical power consists of the torque of each actuators and the rotation angular velocity of the motors. However, the angular velocity range of the joint trajectory planning of the manipulator is relatively small. Therefore, the problem of energy optimization becomes the joint torque optimization problem. For the trajectory optimization procedure, the integral of squared joint torques as the cost function is known from studies in the literature [25]. We introduce the following objective function as a standard for optimization. (9)Econς=∑j=16∫t=0teτj2dt.Here, ς denotes the optimal parameter, which is chosen as the position of mobile platform base and the trajectory of the door angle. To solve the problem of energy consumption optimization during the door-opening procedure, we determine the parameters that minimize Econ:(10)ς∗=argmin⁡Econ.


## 3. Description of the Mobile Manipulator

### 3.1. Kinematics of the Robot

The kinematics model configuration of the 6-DOF Schunk modular manipulator is shown in Figure3.Figure 3
Kinematics configuration of the 6DoF manipulator.The forward kinematics is developed by using the D-H method and represented byT06.(11)T06=R3×3p1×30001=r11r12r13pxr21r22r23pyr31r32r33pz0001,where R3×3 and p1×3 are the rotation and translation matrix. Based on (11), the joint angles by inverse kinematics calculations can be derived as(12)θ1θ2θ3θ4θ5θ6=atan⁡2py,px-atan⁡20,±px2+py2θ23-θ3atan⁡2K,±l32-K2atan⁡2-r13s1+r23c1,-r13c1c23-r23s1c23+r33s23atan⁡2-r13c1c23c4+s1s4-r23s1c23c4-c1s4r33s23c4,-r13c1s23-r23s1s23-r33c23atan⁡2s6,c6,where (13)K=px2+py2+pz2-a22-l322a2,θ23=atan⁡2-a2c3pz-c1px+s1pyl3-a2s3,a2s3-l3pz-a3+a2c3c1px+s1py,ci=cos⁡θi,si=sin⁡θi.Because of the redundancy, there are eight groups of joint angles through the inverse kinematics for giving the position of the end-effector. An algorithm that can minimize the Euler distance in joint space from the initial state should be chosen from the eight solutions [23]. The Jacobin matrix of the 6-DOF manipulator is (14)J=J11J12J13J14J15J16J21J22J23J24J25J26J31J32J33J34J35J36000J44J45J46000J54J550000001.Based on (14), the angular velocity θ˙i and acceleration θ¨i of the joint can be derived as(15)θ˙1θ˙2θ˙3θ˙4θ˙5θ˙6T=J-1P˙θ¨1θ¨2θ¨3θ¨4θ¨5θ¨6T=J-1P¨+J-1P˙,where P is the position and orientation vector of the end-effector of the manipulator.


### 3.2. Dynamics of the Robot

From the inverse kinematics in the previous subsection, the equations for the motion of the system can be written by usingθdt:(16)τ1,τ2,τ3,τ4,τ5,τ6,0T-BF=Ψθd,θ˙d,θ¨d,where F=F1,F2,F3,M1,M2,M3T is the force and torque of the xyz direction used for grasping the door handle. Ψ denotes the dynamics parameter matrix of the door in (16), including the symmetric positive definite manipulator inertia matrix M(θd), the vector of centripetal and Coriolis torques C(θd,θ˙d), and the vector of gravitational torques g(θd). The matrix B can be obtained as (17)BT=∂Φ∂θ,∂Φ∂θd=J,∂Φ∂θd,whereΦ is a function of the relationship between the angle of the door-opening and the joint angle of the manipulator.Φθ,θd=f2(θ)-f1(θd)=0.(18)f1θd=rdcos⁡θdrdsin⁡θdzd,f2θ=xb+d4c2s1s3+c3s1s2+d6s5c1s4-c4s1s2s3-c2c3s1+c5c2s1s3+c3s1s2+a2s1s2yb+d6s5s1s4+c4c1s2s3-c1c2c3-c5c1c2s3+c1c3s2-d4c1c2s3+c1c3s2-a2c1s2zd,where zd is the coordinate value of the z direction in the coordinate system of the door, ci and si are the abbreviations for cos⁡(θi) and sin⁡(θi), respectively, and l1, l2, l3, and l4 are the lengths of the links in Figure 3.We cannot determineτ in (16), which means that the system is indeterminate. F can be measured by using the six-axis F/T sensor, and the forward six rows of (16) can be rewritten as(19)τ-JTF=Mθθ¨+Vθ,θ˙+Gθ,where Mθ is the 6×6 inertia matrix of the manipulator, Vθ,θ˙ is the 6×1 vector of the centrifugal force and Coriolis force, and Gθ is the 6×1 vector of the gravity.


## 3.1. Kinematics of the Robot

The kinematics model configuration of the 6-DOF Schunk modular manipulator is shown in Figure3.Figure 3
Kinematics configuration of the 6DoF manipulator.The forward kinematics is developed by using the D-H method and represented byT06.(11)T06=R3×3p1×30001=r11r12r13pxr21r22r23pyr31r32r33pz0001,where R3×3 and p1×3 are the rotation and translation matrix. Based on (11), the joint angles by inverse kinematics calculations can be derived as(12)θ1θ2θ3θ4θ5θ6=atan⁡2py,px-atan⁡20,±px2+py2θ23-θ3atan⁡2K,±l32-K2atan⁡2-r13s1+r23c1,-r13c1c23-r23s1c23+r33s23atan⁡2-r13c1c23c4+s1s4-r23s1c23c4-c1s4r33s23c4,-r13c1s23-r23s1s23-r33c23atan⁡2s6,c6,where (13)K=px2+py2+pz2-a22-l322a2,θ23=atan⁡2-a2c3pz-c1px+s1pyl3-a2s3,a2s3-l3pz-a3+a2c3c1px+s1py,ci=cos⁡θi,si=sin⁡θi.Because of the redundancy, there are eight groups of joint angles through the inverse kinematics for giving the position of the end-effector. An algorithm that can minimize the Euler distance in joint space from the initial state should be chosen from the eight solutions [23]. The Jacobin matrix of the 6-DOF manipulator is (14)J=J11J12J13J14J15J16J21J22J23J24J25J26J31J32J33J34J35J36000J44J45J46000J54J550000001.Based on (14), the angular velocity θ˙i and acceleration θ¨i of the joint can be derived as(15)θ˙1θ˙2θ˙3θ˙4θ˙5θ˙6T=J-1P˙θ¨1θ¨2θ¨3θ¨4θ¨5θ¨6T=J-1P¨+J-1P˙,where P is the position and orientation vector of the end-effector of the manipulator.


## 3.2. Dynamics of the Robot

From the inverse kinematics in the previous subsection, the equations for the motion of the system can be written by usingθdt:(16)τ1,τ2,τ3,τ4,τ5,τ6,0T-BF=Ψθd,θ˙d,θ¨d,where F=F1,F2,F3,M1,M2,M3T is the force and torque of the xyz direction used for grasping the door handle. Ψ denotes the dynamics parameter matrix of the door in (16), including the symmetric positive definite manipulator inertia matrix M(θd), the vector of centripetal and Coriolis torques C(θd,θ˙d), and the vector of gravitational torques g(θd). The matrix B can be obtained as (17)BT=∂Φ∂θ,∂Φ∂θd=J,∂Φ∂θd,whereΦ is a function of the relationship between the angle of the door-opening and the joint angle of the manipulator.Φθ,θd=f2(θ)-f1(θd)=0.(18)f1θd=rdcos⁡θdrdsin⁡θdzd,f2θ=xb+d4c2s1s3+c3s1s2+d6s5c1s4-c4s1s2s3-c2c3s1+c5c2s1s3+c3s1s2+a2s1s2yb+d6s5s1s4+c4c1s2s3-c1c2c3-c5c1c2s3+c1c3s2-d4c1c2s3+c1c3s2-a2c1s2zd,where zd is the coordinate value of the z direction in the coordinate system of the door, ci and si are the abbreviations for cos⁡(θi) and sin⁡(θi), respectively, and l1, l2, l3, and l4 are the lengths of the links in Figure 3.We cannot determineτ in (16), which means that the system is indeterminate. F can be measured by using the six-axis F/T sensor, and the forward six rows of (16) can be rewritten as(19)τ-JTF=Mθθ¨+Vθ,θ˙+Gθ,where Mθ is the 6×6 inertia matrix of the manipulator, Vθ,θ˙ is the 6×1 vector of the centrifugal force and Coriolis force, and Gθ is the 6×1 vector of the gravity.


## 4. Optimization Method

In this part, we focus on the optimal method for the door-opening procedure. Sinceθdt is an infinite dimensional parameter, it is difficult to find the optimal solution that minimizes Econ rigorously. Therefore, we approximate θdt as the fifth-order spline functions of time and find the coefficients of splines that minimize the cost function. The position of the mobile manipulators base (xb,yb) is discretized into a grid, and at each grid point, the quasi-optimal motion of the door is calculated by using the spline functions.


### 4.1. Optimal Motion of Door

We divide the time interval[0,te] by n and assume that the trajectory of θdt in each time interval [ti,ti+1]  (i=0,…,n-1 and tj=jte/n for j=0,…,n) is expressed by a fifth-order polynomial function of time, Xit, as(20)Xit=θdi+a1it-ti+a2it-ti2+a3it-ti3+a4it-ti4+a5it-ti5,where a1i, a2i, a3i, a4i, and a5i are the coefficients of the polynomial. In order to make the input torque τ continuous, we choose the function such that Xi(ti+1)=Xi+1(ti+1),  X˙i(ti+1)=X˙i+1(ti+1),  X¨i(ti+1)=X¨i+1(ti+1), for i=0,…,n-2, and (20) satisfies (2). The coefficients a1i, a2i, and a3i are constant. Therefore, there are 3n-1 independent parameters, and they can be shown as(21)Ζ=θd1,…,θdn-1,a40,…,a4n-1,a50,…,a5n-1.We assume that the door angleθdt is a monotonically increasing function at t=ti with constraint 0≤θd0≤⋯≤θdn-1≤θdte. At each position of the mobile manipulator’s base, we search for the values of Ζ that minimize the Econ by using the Quasi-Newton method. We can find a vector h=h1,h2,h3,h4,h5,h6,h7T that conforms to hTB=0. Eq. (16) is multiplied by the vector h, and thus, we obtain(22)h1τ1+h2τ2+h3τ3+h4τ4+h5τ5+h6τ6=hTΨ.The manipulator joint torqueτ satisfying (22) can obtain a specified motion of the door θd(t). When τ is minimized, τ can be shown as(23)τ=kτh1,h2,h3,h4,h5,h6T,where kτ is a scalar parameter. By using (23), we can calculate τ uniquely. We propose a method, as shown in Algorithm 1, to acquire the optimal motion of the door and the corresponding torque of each joint.Algorithm 1:Select the optimal motion of door.
Input: the spline functions Xi(t) for all time internals Output: the torques τ and the parameters of Eq. (21) For i in range (1,n) Based on Quasi-Newton method with constraintmin⁡Econ ObtainZi=(θdi,a4i,a5i) solveθ1i,θ2i,…,θ6i=InverseKinematicsXitτ1i,τ2i,…,τ6i=eq.23θ1i,θ2i,…,θ6i,θdi τ=sum(τ(i)) i=i+1 End Return τ, Zi=(θdi,a4i,a5i)


### 4.2. Optimal Position of the Mobile Manipulator’s Base

We determine the optimal position of the mobile manipulator’s base and the end-effector’s grasp by using the exhaustive method. The region of(xb,yb), defined by xbmin,xbmax×ybmin,ybmax, is divided into a grid, where each rectangle is given by Δx×Δy. By calculating the objective function at each grid point using the method in Section 4.1, with a set of [τi,Xbi,Ybi], we can determine the optimal position of the mobile manipulator’s base. The scheme for optimizing the position is shown in Algorithm 2.Algorithm 2:Search the optimal position of robot base.
Input: [ τ i , X b i , Y b i ] Output:the minimum Econ and the corresponding base position [Xbmin,Ybmin] Econmin=∞ For i in range (1,n) Econi IfEconi<Econmin thenEconmin=Econi,  Xbmin=Xbi,  Ybmin=Ybi, i=i+1 End Return E c o n  m i n ,   ( X b m i n , Y b m i n )


### 4.3. Numerical Simulations

In this section, we describe the acquisition of the optimal solution by using numerical simulations. We used MATLAB to find the optimal values ofZ under the constraint that 0≤θd0≤⋯≤θdn-1≤θdt. The open angle of the door, θd(t), was set to be π/3 at te=20 s. The time of the door pulling interval [0,20] was divided into ten subintervals; that is, n=10. The measured length of the door was 0.85 [m]. The mass and inertia of the door were set to 35.7 [kg] and 8.598 [kg·m2], respectively. The lengths of l1, l2, l3, and l4 shown in Figure 3 were 0.3 [m], 0.3 [m], 0.305 [m], and 0.415 [m], respectively. The mass and inertia of the manipulator were as provided in [26]. In the simulation calculation, we chose the grid points of (xb,yb), that is, Δx and Δy as 0.02 [m], the search area of xbmin,xbmax×ybmin,ybmax was 0.7m,1.1m×0.3m,0.9m, and the objective function defined by (9) was calculated at each point. Figure 4 shows the contour plot of the objective function Econ. It shows that the position of the minimum Econ is (xb,yb)=(0.82m,0.48m) and the value of Econ is 8.28×105 [N2m2s]. Figure 5 shows the motion process of the door in the planning time.Figure 4
Objective for base position(xb,yb).Figure 5
Optimal moment of the door.


## 4.1. Optimal Motion of Door

We divide the time interval[0,te] by n and assume that the trajectory of θdt in each time interval [ti,ti+1]  (i=0,…,n-1 and tj=jte/n for j=0,…,n) is expressed by a fifth-order polynomial function of time, Xit, as(20)Xit=θdi+a1it-ti+a2it-ti2+a3it-ti3+a4it-ti4+a5it-ti5,where a1i, a2i, a3i, a4i, and a5i are the coefficients of the polynomial. In order to make the input torque τ continuous, we choose the function such that Xi(ti+1)=Xi+1(ti+1),  X˙i(ti+1)=X˙i+1(ti+1),  X¨i(ti+1)=X¨i+1(ti+1), for i=0,…,n-2, and (20) satisfies (2). The coefficients a1i, a2i, and a3i are constant. Therefore, there are 3n-1 independent parameters, and they can be shown as(21)Ζ=θd1,…,θdn-1,a40,…,a4n-1,a50,…,a5n-1.We assume that the door angleθdt is a monotonically increasing function at t=ti with constraint 0≤θd0≤⋯≤θdn-1≤θdte. At each position of the mobile manipulator’s base, we search for the values of Ζ that minimize the Econ by using the Quasi-Newton method. We can find a vector h=h1,h2,h3,h4,h5,h6,h7T that conforms to hTB=0. Eq. (16) is multiplied by the vector h, and thus, we obtain(22)h1τ1+h2τ2+h3τ3+h4τ4+h5τ5+h6τ6=hTΨ.The manipulator joint torqueτ satisfying (22) can obtain a specified motion of the door θd(t). When τ is minimized, τ can be shown as(23)τ=kτh1,h2,h3,h4,h5,h6T,where kτ is a scalar parameter. By using (23), we can calculate τ uniquely. We propose a method, as shown in Algorithm 1, to acquire the optimal motion of the door and the corresponding torque of each joint.Algorithm 1:Select the optimal motion of door.
Input: the spline functions Xi(t) for all time internals Output: the torques τ and the parameters of Eq. (21) For i in range (1,n) Based on Quasi-Newton method with constraintmin⁡Econ ObtainZi=(θdi,a4i,a5i) solveθ1i,θ2i,…,θ6i=InverseKinematicsXitτ1i,τ2i,…,τ6i=eq.23θ1i,θ2i,…,θ6i,θdi τ=sum(τ(i)) i=i+1 End Return τ, Zi=(θdi,a4i,a5i)


## 4.2. Optimal Position of the Mobile Manipulator’s Base

We determine the optimal position of the mobile manipulator’s base and the end-effector’s grasp by using the exhaustive method. The region of(xb,yb), defined by xbmin,xbmax×ybmin,ybmax, is divided into a grid, where each rectangle is given by Δx×Δy. By calculating the objective function at each grid point using the method in Section 4.1, with a set of [τi,Xbi,Ybi], we can determine the optimal position of the mobile manipulator’s base. The scheme for optimizing the position is shown in Algorithm 2.Algorithm 2:Search the optimal position of robot base.
Input: [ τ i , X b i , Y b i ] Output:the minimum Econ and the corresponding base position [Xbmin,Ybmin] Econmin=∞ For i in range (1,n) Econi IfEconi<Econmin thenEconmin=Econi,  Xbmin=Xbi,  Ybmin=Ybi, i=i+1 End Return E c o n  m i n ,   ( X b m i n , Y b m i n )


## 4.3. Numerical Simulations

In this section, we describe the acquisition of the optimal solution by using numerical simulations. We used MATLAB to find the optimal values ofZ under the constraint that 0≤θd0≤⋯≤θdn-1≤θdt. The open angle of the door, θd(t), was set to be π/3 at te=20 s. The time of the door pulling interval [0,20] was divided into ten subintervals; that is, n=10. The measured length of the door was 0.85 [m]. The mass and inertia of the door were set to 35.7 [kg] and 8.598 [kg·m2], respectively. The lengths of l1, l2, l3, and l4 shown in Figure 3 were 0.3 [m], 0.3 [m], 0.305 [m], and 0.415 [m], respectively. The mass and inertia of the manipulator were as provided in [26]. In the simulation calculation, we chose the grid points of (xb,yb), that is, Δx and Δy as 0.02 [m], the search area of xbmin,xbmax×ybmin,ybmax was 0.7m,1.1m×0.3m,0.9m, and the objective function defined by (9) was calculated at each point. Figure 4 shows the contour plot of the objective function Econ. It shows that the position of the minimum Econ is (xb,yb)=(0.82m,0.48m) and the value of Econ is 8.28×105 [N2m2s]. Figure 5 shows the motion process of the door in the planning time.Figure 4
Objective for base position(xb,yb).Figure 5
Optimal moment of the door.


## 5. Experimental Results

### 5.1. Robot System

The mobile robot system (see Figure6) consists of a 6-DOF modular manipulator produced by the SCHUNK Company, a four-wheeled platform with two drive wheels, and a two free wheels and two-finger gripper mounted on the arm wrist module. The 6-DOF modular manipulator comprises three types of modules, and each joint module consists of a brushless DC motor, a harmonic drive, a braking system, and an encoder. The mobile robot system is equipped with various types of sensors, including two six-axis force sensors and a camera. The two six-axis force sensors are mounted on the manipulator wrist module and the base module, respectively, to measure the mechanics date of the door-opening procedure. A joystick with associated force feedback control from the base-mounted 6-axis force/torque sensor is used for teleoperation. A camera is mounted on top of the frame to support 3D displays for continuous teleportation and object recognition. An API T3 (Automated Precision Inc.) laser tracker system was used to measure the base position in the experiment.Figure 6
Software and hardware of the mobile manipulator.


### 5.2. Controller Description

The structure of the proposed controller scheme is illustrated in Figure7. The control system is shown in terms of the composition of the PD-computed torque controller. The computed torque controller (CTC) uses a feedback linearization method. It is assumed that the desired motion trajectory for the manipulator is determined by a path planner [27]. The tracking error is defined as (24)et=θrt-θt,where e(t) is the error of the plant, θr(t) is the desired input variable, which is the desired joint angular displacement in the control system, and θ(t) is the actual joint angular displacement.Figure 7
Block diagram of the computed torque controller.The control torque is described as(25)τ~=Vθ,θ˙θ˙+Gθ+Mθu.The dynamic model in (19) is equivalent to a decoupled linear time-invariant system:(26)θ¨=u.Considering that the desired trajectoryθr(t) is determined, θ˙r and θ¨r are known. This is a nonlinear feedback control law that guarantees tracking of the desired door-opening trajectory. Selecting PD feedback for u(t) results in the PD-computed torque controller.(27)u=θ¨r+Kdθ˙r-θ˙+Kpθd-θ=θ¨r+Kde˙+Kpe,where Kd and Kp are the controller gains and are positive definite diagonal matrices, Kp=diag(8.2,10.4,8.5,7.4,12.5,9.0), Kd=diag(7.0,9.5,8.0,7.0,11.5,8.5). The closed-loop system equation is(28)θ¨r+Kde˙+Kpe=0.By substituting (27) into (25), we obtain the complete expression of the control law:(29)τ~=Mθθ¨r+Kde˙+Kpe+Vθ,θ˙θ˙+Gθ.


### 5.3. Experimental Setup

The experiments were conducted in a simulated NPPs internal environment and a real NPPs fire door with a door closer was used. The trajectory during door-opening was calculated using the method described in Section2.2. The entire door-opening procedure was realized using C++ programming with a cycle of 20 ms. Four groups of experiments using different base positions were conducted. In the first group, the base position, denoted by G1, was (xb=0.82 m, yb=0.48 m). This is the point of the minimum objective function according to the numerical simulations. In the second group, denoted by G2, was (xb=0.85 m, yb=0.55 m). In the third group, the base position, denoted by G3, was (xb=0.72 m, yb=0.55 m). In the fourth group, the base position, denoted by G4, was (xb=0.95 m, yb=0.7 m). The points G1, G2, G3, and G4 have been pointed and marked in Figure 4. G1 is the point of the minimum objective function as shown in Figure 4.Figure8 shows sequential pictures of the door-opening experiments. The duration of the door-opening procedure is 20 s, which is the same as that the numerical simulation described in Section 4.3. The mobile robot system can measure the instantaneous current in the door-opening procedure. The torque of each joint was calculated using the values of the instantaneous current and (8). Therefore, the value of Econ could be calculated in the various groups of experiments. The test data were acquired at a sampling frequency of 50 Hz. The experiments were conducted more than three times under the same conditions to ensure the repeatability and effectiveness of the experimental results. During the comparison of the energy consumption in the different groups of experiments, the end-effector grasped the door handle in the same positions according to the door-opening path planning provided in Section 2.Figure 8
Door-opening experiment: (a) initial moments; (b) pulling the door.
(a) (b)The performance parameters of the motors and actuating devices of joints are shown in Table1.Table 1
Performance parameters of drive devices and actuating devices of joints.
Joint of armj U i n/(V) K m j/(mN⋅m) i j R j / ( Ω ) η j / ( % ) Joints 1 and 2 24 31.4 596 1.5 0.85 Joints 3 and 4 24 38 625 0.6 0.85 Joints 5 and 6 24 16 552 0.1 0.85


### 5.4. Comparison of the Different Base Positions

The matrixτj when the robot base is positioned at G1, G2, G3, and G4 can be calculated by using (8) and the performance parameters of the joints, shown in Table 1. The values of ∑j=16τj2 during the door-opening procedure are shown in Figure 9. The sums of the objective function for the different base positions are shown in Table 2.Table 2
Value ofEcon at different base positions.
Objective function G1 G2 G3 G4 E c o n   ( N 2 m 2 s ) 2.97 × 1 0 6 1.65 × 1 0 7 2.79 × 1 0 7 2.95 × 1 0 7Figure 9
Comparison of the square sums of joint torque at the different base positions.TheEcon curves shown in Figure 9 demonstrate that energy consumption is optimized at base position G1. Figure 9 shows that the torque is increased during the start and end of the door-opening procedure because of the resistance of the door handle and the door closer. It can be clearly seen in Figure 9 that, in comparison with other selected positions, over time G1 is better than the other selected base positions minimizing the objective function, not only in terms of the accumulated value of Econ, shown in Table 2, but also in terms of sustaining the lowest value of ∑j=16τj2. This verifies the effectiveness of the proposed algorithm. The experiment value of Econ is larger than the simulation value at the same base position (G1), because the simulation only considered the mechanical power of (6); however, the experiment value included the generated mechanical power, heat power, and power losses due to factors such as friction.Figures10 and 11 show the optimal trajectories of the joint angle θi and the torque τ of each joint for G1. It can be seen in Figure 11 that the values of the joints’ torques τ2,τ3 are larger than those of the other torques. Figure 10 shows that in the motion during the door-opening procedure the angle degree of Joints 2 and 3 is larger than that of the other joints.Figure 10
Angle of each joint (G1).Figure 11
Torque of each joint (G1).


## 5.1. Robot System

The mobile robot system (see Figure6) consists of a 6-DOF modular manipulator produced by the SCHUNK Company, a four-wheeled platform with two drive wheels, and a two free wheels and two-finger gripper mounted on the arm wrist module. The 6-DOF modular manipulator comprises three types of modules, and each joint module consists of a brushless DC motor, a harmonic drive, a braking system, and an encoder. The mobile robot system is equipped with various types of sensors, including two six-axis force sensors and a camera. The two six-axis force sensors are mounted on the manipulator wrist module and the base module, respectively, to measure the mechanics date of the door-opening procedure. A joystick with associated force feedback control from the base-mounted 6-axis force/torque sensor is used for teleoperation. A camera is mounted on top of the frame to support 3D displays for continuous teleportation and object recognition. An API T3 (Automated Precision Inc.) laser tracker system was used to measure the base position in the experiment.Figure 6
Software and hardware of the mobile manipulator.


## 5.2. Controller Description

The structure of the proposed controller scheme is illustrated in Figure7. The control system is shown in terms of the composition of the PD-computed torque controller. The computed torque controller (CTC) uses a feedback linearization method. It is assumed that the desired motion trajectory for the manipulator is determined by a path planner [27]. The tracking error is defined as (24)et=θrt-θt,where e(t) is the error of the plant, θr(t) is the desired input variable, which is the desired joint angular displacement in the control system, and θ(t) is the actual joint angular displacement.Figure 7
Block diagram of the computed torque controller.The control torque is described as(25)τ~=Vθ,θ˙θ˙+Gθ+Mθu.The dynamic model in (19) is equivalent to a decoupled linear time-invariant system:(26)θ¨=u.Considering that the desired trajectoryθr(t) is determined, θ˙r and θ¨r are known. This is a nonlinear feedback control law that guarantees tracking of the desired door-opening trajectory. Selecting PD feedback for u(t) results in the PD-computed torque controller.(27)u=θ¨r+Kdθ˙r-θ˙+Kpθd-θ=θ¨r+Kde˙+Kpe,where Kd and Kp are the controller gains and are positive definite diagonal matrices, Kp=diag(8.2,10.4,8.5,7.4,12.5,9.0), Kd=diag(7.0,9.5,8.0,7.0,11.5,8.5). The closed-loop system equation is(28)θ¨r+Kde˙+Kpe=0.By substituting (27) into (25), we obtain the complete expression of the control law:(29)τ~=Mθθ¨r+Kde˙+Kpe+Vθ,θ˙θ˙+Gθ.


## 5.3. Experimental Setup

The experiments were conducted in a simulated NPPs internal environment and a real NPPs fire door with a door closer was used. The trajectory during door-opening was calculated using the method described in Section2.2. The entire door-opening procedure was realized using C++ programming with a cycle of 20 ms. Four groups of experiments using different base positions were conducted. In the first group, the base position, denoted by G1, was (xb=0.82 m, yb=0.48 m). This is the point of the minimum objective function according to the numerical simulations. In the second group, denoted by G2, was (xb=0.85 m, yb=0.55 m). In the third group, the base position, denoted by G3, was (xb=0.72 m, yb=0.55 m). In the fourth group, the base position, denoted by G4, was (xb=0.95 m, yb=0.7 m). The points G1, G2, G3, and G4 have been pointed and marked in Figure 4. G1 is the point of the minimum objective function as shown in Figure 4.Figure8 shows sequential pictures of the door-opening experiments. The duration of the door-opening procedure is 20 s, which is the same as that the numerical simulation described in Section 4.3. The mobile robot system can measure the instantaneous current in the door-opening procedure. The torque of each joint was calculated using the values of the instantaneous current and (8). Therefore, the value of Econ could be calculated in the various groups of experiments. The test data were acquired at a sampling frequency of 50 Hz. The experiments were conducted more than three times under the same conditions to ensure the repeatability and effectiveness of the experimental results. During the comparison of the energy consumption in the different groups of experiments, the end-effector grasped the door handle in the same positions according to the door-opening path planning provided in Section 2.Figure 8
Door-opening experiment: (a) initial moments; (b) pulling the door.
(a) (b)The performance parameters of the motors and actuating devices of joints are shown in Table1.Table 1
Performance parameters of drive devices and actuating devices of joints.
Joint of armj U i n/(V) K m j/(mN⋅m) i j R j / ( Ω ) η j / ( % ) Joints 1 and 2 24 31.4 596 1.5 0.85 Joints 3 and 4 24 38 625 0.6 0.85 Joints 5 and 6 24 16 552 0.1 0.85


## 5.4. Comparison of the Different Base Positions

The matrixτj when the robot base is positioned at G1, G2, G3, and G4 can be calculated by using (8) and the performance parameters of the joints, shown in Table 1. The values of ∑j=16τj2 during the door-opening procedure are shown in Figure 9. The sums of the objective function for the different base positions are shown in Table 2.Table 2
Value ofEcon at different base positions.
Objective function G1 G2 G3 G4 E c o n   ( N 2 m 2 s ) 2.97 × 1 0 6 1.65 × 1 0 7 2.79 × 1 0 7 2.95 × 1 0 7Figure 9
Comparison of the square sums of joint torque at the different base positions.TheEcon curves shown in Figure 9 demonstrate that energy consumption is optimized at base position G1. Figure 9 shows that the torque is increased during the start and end of the door-opening procedure because of the resistance of the door handle and the door closer. It can be clearly seen in Figure 9 that, in comparison with other selected positions, over time G1 is better than the other selected base positions minimizing the objective function, not only in terms of the accumulated value of Econ, shown in Table 2, but also in terms of sustaining the lowest value of ∑j=16τj2. This verifies the effectiveness of the proposed algorithm. The experiment value of Econ is larger than the simulation value at the same base position (G1), because the simulation only considered the mechanical power of (6); however, the experiment value included the generated mechanical power, heat power, and power losses due to factors such as friction.Figures10 and 11 show the optimal trajectories of the joint angle θi and the torque τ of each joint for G1. It can be seen in Figure 11 that the values of the joints’ torques τ2,τ3 are larger than those of the other torques. Figure 10 shows that in the motion during the door-opening procedure the angle degree of Joints 2 and 3 is larger than that of the other joints.Figure 10
Angle of each joint (G1).Figure 11
Torque of each joint (G1).


## 6. Conclusion

In this paper, we proposed a novel energy consumption optimization scheme for a mobile manipulator executing the door-opening task. Focusing on the power consumption of the manipulator during the entire task period, we chose the quadric terms of the joint torques as the objective function. Furthermore, in this study a two-step optimization procedure was developed to solve the corresponding joint trajectories of the manipulator. In the first step, the feasible base positions of the manipulator are decentralized into a grid in order to simplify the entire optimization process. A piecewise fifth-order polynomials over time is utilized to parameterize the rotational motion of the door. By applying the Quasi-Newton method, the local optimal trajectories of the manipulator are obtained for a given base position. In the second step, the optimal base position is attained via searching the decentralized grid of the feasible base positions. The numerical results when the proposed method was applied showed that energy consumption was optimized at base position G1. The experimental results for door-opening at the different base positions demonstrate the effectiveness of the method proposed in this study. The proposed method will be useful for the development of the NPP rescue robots.



---

*Source: 8987953-2018-06-25.xml*