Bilateral teleoperation of stochastic port‐Hamiltonian systems using energy tanks

In this article we consider the general problem of how to properly endow a stochastic port‐Hamiltonian system (SPHS) with an energy tank, that is an energy reservoir that allows to guarantee the passivity property. We show that a stochastic bilateral teleoperation system, composed by a master robot and a slave robot modeled as SPHS, can be connected in a power‐preserving manner to energy tanks. The stored energy is continuously monitored to keep the system passive despite time‐varying communication delays and interaction with unknown environment that may destabilize the overall system. We will address latter problem considering a SPHS affected by a noise composed by a linear, multiplicative component in Itô form plus an additive one. We underline that such a scenario requires the introduction of an ad hoc notion of passivity.


INTRODUCTION
In recent decades there has seen a constant growing interest in bilateral teleoperation. This scenario is characterized by a user acting on a haptic device, usually referred to as master robot, to remotely control a slave robot. Typical examples of such systems can be found in robotic minimally-invasive surgery R-MIS, [1][2][3] or concerning maintenance in hazard and harsh environments. 4,5 The proper control of a teleoperated system is not a trivial task. In fact, from a purely theoretical control point of view, the master and the slave should be connected in such a way as to ensure the overall stability despite several factors concurring to potentially destabilize the whole system. 6 Among such sources of instability a relevant role is played by the unknown environment the slave robot is interacting with, as well as by the time delay in the communication channel connecting the master and the slave side. For instance in Reference 7 a time varying stochastic delay has been considered.
Within the teleoperation research topic, the Passivity theory approach has proved to be particularly suitable for addressing the stability problem of teleoperated systems. 6 Since, in general, it is not possible to mathematically describe the environment the slave robot will interact with and the operator, it is common to require the passivity assumption on such systems. Then, control strategies are derived to ensure stability of the whole system. 6,8 Moreover, when dealing with real world applications, the communication channel connecting the master and slave robot is subjected to unknown and time-varying delays. As for the case of unknown exogenous systems, passivity theory can be exploited to derive conditions for designing control architectures that guarantee a stable behavior.
Recently, a new approach has been introduced to enforce passivity by endowing the teleoperation system with two energy tanks/reservoirs (at the master and at the slave side) where to store the energy dissipated by the system. This energy can be used to modulate in a passive way commands that would have otherwise led to unstable behaviors for the system,. 3,5,9,10 Such an approach can deal with time-varying delays providing less conservative control strategies than the ones developed using, for example, adaptive proportional-derivative (PD) based controls 11 or wave variables. 12 It is worth mentioning that the energy tank approach is strictly related to the mathematical theory of port-Hamiltonian systems (PHS), a well established theory that merges geometric mechanics with bond graph theory, focusing on the overall energy of the system. 13 In particular, when considering teleoperated systems and energy tanks, PHS allow to study the resulting system in an abstract framework, proving general results both on the passivity and the stability of the system.
The main contribution of the present work is to formally extend the energy tank approach by considering proper stochastic perturbations of PHS. To the best of our knowledge, very few results have been already obtained about stochastic PHS (SPHS). [14][15][16][17] When a noise enters the system, passivity becomes a more subtle and delicate point than within the standard deterministic setting. 14 In fact, the noise can inject energy into the system destroying passivity. In the present work we consider the input-state-output PHS framework showing how an energy tank, that stores energy dissipated by the system plus expected energy coming from the noise term, should be added to the system to passify the control strategy.
In the first part of the present article we show how to derive a formal stochastic energy tank for a SPHS with multiplicative linear noise of Brownian type. Since, almost surely, the Brownian motion has not absolutely continuous paths, the standard notion of integration cannot be used, hence an ad hoc concept has to be introduced. To this end, over the years, two notions have been proven to be particularly effective: the Itô stochastic integral and the Stratonovich one. 18 Such notions of integrals lead to different results, because of the different way they exploit the flux of information characterizing the studied system. Concerning our analysis, we will consider the stochastic integration in the Itô sense, showing further how energy tank and passivity results change when Stratonovich calculus is used instead. For the sake of completeness, we have recalled the main differences between these two types of stochastic integrals in Appendix A.
In the second part of our article, we consider the case of a SPHS with additive noise. From a stability point of view, an additive noise poses several difficulties. In particular, when a noise is added to a deterministic system with stable solution X * and we still aim at maintaining the system stability, we have to consider a noise that vanishes at the stable solution X * . Otherwise the noise will constantly force the system away from the stable point. Since passivity is strongly related to stability, the instability generated by a non-vanishing noise affects also the passivity of the system. It follows that one of the most difficult case to treat is the one characterized by additive noise since it does not vanish at the stable solution. To overcome this problem, more general types of stability should be introduced. 16 A process is called stable if it possesses a long-time stable distribution. In the present work, we consider a weak notion of passivity, called ultimately passivity, to show how the passivity problem can be tackled in the presence of an additive noise. Such a notion is related to the one of ultimately bounded stability. 19 We also recall that the mentioned weaker notion of passivity has already been introduced, under the name of weak stochastic passivity, in Reference 20. We decided to chose a different name to better highlight the connection between the notions of ultimately passive and of ultimately bounded stability, this last concept aiming at studying the ergodic behavior of the system. We will exploit the ultimately bounded stability definition to prove that a stochastic system characterized by an additive noise is ultimately passive when endowed with an energy tank.
In the third part of the article, we use the results obtained for SPHS endowed with tanks and affected by multiplicative/additive stochastic noise in the context of bilateral teleoperation with communication delay. The case of the additive noise appears of particular interest. From a mathematical point of view, previous setting translates into studying and characterizing the transition semigroup of delayed systems with degenerate noise. Very few results exist that can handle both delay and degenerate noise, moreover a systematic treatment is completely missing. A general and detailed treatment on Feller property and related regularity characteristics of the transition semigroup can be found in References 21-26. Moreover, in Reference 26, the authors consider the case of systems of oscillators, allowing the noise to be degenerate. By taking into account stochastic noises within the tank-based two-layer approach, 5 we will prove weak passivity as well as existence and uniqueness of an invariant measure and convergence of the system towards the (unique) invariant distribution. We refer the interested reader to References 21 and 22 for an in-depth treatment of ergodicity and Feller property for infinite dimensional and delayed stochastic systems, or to Reference 27 for an extensive overview of stability for stochastic delay equations.
For the sake of completeness, let us remark that explicit derivations of invariant measure for stochastic systems are rare and difficult to obtain. In fact, besides the class of systems usually referred to as gradient systems, it is not possible to obtain the invariant measure in closed form. Within the framework considered in the present work, latter task becomes even more difficult. Indeed, we are dealing with a time-delaying and degenerate system. Nonetheless, several results regarding absolute continuity, regularity and integrability properties, already exist in literature. Therefore, being beyond the scope of the present research, we will not enter into details regarding the explicit form of the limiting invariant measure, rather limiting ourselves to prove existence and uniqueness for an invariant measure characterizing the system we are interested in.
The innovative contributions of the present article are: • The first formal treatment of stochastic PHS endowed with energy tank affected by multiplicative and additive noise; • Existence and uniqueness of an invariant measure for a delayed stochastic system with degenerate noise; • Analysis of stochastic passivity and stability for a bilateral teleoperation system with both stochastic noise and stochastic communication delay.
The article is structured as follows: in Section 2 we introduce the general setting of stochastic PHS according to the Itô setting, considering both multiplicative and additive noise. Moreover a suitable weak notion of stochastic passivity is also introduced. In Section 3 a stochastic bilateral teleoperated system with communication delays is described and the two-layer control architecture is re-written in the proposed stochastic scenario. Sections 4 and 5 are focused on the case of multiplicative noise and of additive noise, respectively. Then, Section 6 presents a numerical example, whereas conclusions and perspectives are drawn in Section 7. Appendix A recalls both the Itô and the Stratonovich approaches to stochastic integration. It is worth highlighting that in the additive noise case the two types of stochastic integrals do coincide.

STOCHASTIC PORT-HAMILTONIAN SYSTEMS AND ENERGY TANKS
The goal of this section is to show how properly endowing a stochastic PHS with an energy tank. Details for PHS, such as passivity analysis and energy tank in the deterministic case, can be found in many references as, for example, References 3,13. As mentioned above, when a stochastic perturbation affects a deterministic PHS, passivity analysis becomes a more delicate point to address. Let Σ be a Stochastic Port-Hamiltonian System (SPHS) in input-state-output form 13 perturbed by a general noise of Brownian type 14,15 where X ∈ R n , J = −J T and R ≻ 0 are two given matrices of suitable dimensions; H ∶ R n → R is the Hamiltonian function, representing the total energy of the system; u ∈  represents the signal for the control input, while y ∈  ∶=  * is the output (as usual  * denotes the dual of the space  ); x denotes the gradient operator with respect to x; ∈ R n×d , with W(t) being a d-dimensional Brownian motion adapted to a reference filtered probability space ( Ω,  , ( t ) t≥0 , P ) satisfying usual assumptions, namely right-continuity and saturation by P-null sets; dW(t) recalls us that we are considering the integration in the sense of Itô. 18 In subsequent sections, we shall also consider the Stratonovich approach to stochastic integration.
As usual in stochastic analysis, we will use the convention that X denotes a random variable where x represents its expected (deterministic) value conditional to  t , representing the information available up to time t. We now introduce the definitions of stochastic passivity and the weaker notion of ultimately stochastic passivity for stochastic systems, 15 that will be proved for SPHS Σ with multiplicative and additive noise, respectively. Definition 1 (Stochastic passivity). The non-linear controlled stochastic system of the form { dX(t) = (X(t), u(t))dt + (X(t))dW(t) is stochastic passive with respect to a storage function V(t, x) if for any (x, t) it holds where  is the infinitesimal generator associated to the system (2) where we have denoted by Tr the trace of a matrix.
Definition 2 (Ultimately stochastic passivity and Strictly ultimately stochastic passivity). The non-linear controlled stochastic system of the form { dX(t) = (X(t), u(t))dt + (X(t))dW(t), is said to be ultimately stochastic passive with respect to a storage function V(t, x) if for x C ∈ R n and for any (x, t) such that ||x − x C || ≥ C, for a given constant C > 0 called passivity radius, it holds being  the infinitesimal generator If further there exists C > 0 such that for ||x − x C || ≥ C it holds then system (3) is said to be strictly ultimately stochastic passive.
Remark 1. We adopted the "new" name of ultimately stochastic passivity, even though the related definition coincides with the one given in Reference 20, to highlight its relationship with the more general notion of ultimately bounded convergence in the deterministic setting. 19 In fact, in the presence of a non-vanishing term, as in the present context with an additive noise, the process does not converge to an equilibrium solution. In particular, in the deterministic setting, the bounded stability can be proven to hold if there exists a Lyapunov function V , so thatV(x) < 0, ∀ x such that ||x − x C || > C, with a given x C . In the stochastic setting we could retrieve a similar result choosing as candidate Lyapunov function the Hamiltonian of the system. If the process X is ultimately stochastic passive, then for the autonomous process with null control, that is, u ≡ 0, it holds V(x) < − < 0, ∀ x such that ||x − x C || > C.

Energy tank for stochastic port-Hamiltonian systems with multiplicative noise
The SPHS Σ m has a multiplicative noise if (x) is a n × d−dimensional matrix which depends on x with the property that, given a stationary solution X * , it must hold (X * ) = 0. Then, we can prove the following proposition.

Proposition 1. The SPHS Σ m in (4) is stochastic passive if and only if
Proof. The infinitesimal generator  of the Itô process (1), is therefore if condition (5) holds we have which is the definition of stochastic passivity. ▪ Assuming that the SPHS with multiplicative noise (4) is stochastic passive, then we can endow the system with a (virtual) energy tank which allows to store the energy dissipated by Σ m and to manage the energy exchanged at its power port with other subsystems. We will see in the next section how tanks are useful to guarantee passivity (and so stability) for teleoperated systems by modulating the control commands as already done for deterministic PHS. 3 Let (x) > 0 be the mean power dissipated by the stochastic system (4) where D D (x) is the contribution of the deterministic part, while D S (x) represents the stochastic contribution. The mathematical model of the tank storing the dissipated energy (6) is where x E is the state variable of the tank and the energy within the tank at time t is computed as The dynamical equation in (7) can be derived by observing that The parameter ∈ {0, 1} is used to able or disable the storing of the dissipated energy in case the energy in the tank reaches an application-dependent upper bound T max set for safety reason. 3 The tank has to be initialized therefore T(x E (0)) > T min and energy extraction is prohibited if T(x E (t)) < T min . The energetic port of the tank (u E , y E ) is connected to the energetic port (u, y) of the SPHS Σ m via a power preserving interconnection that should guarantee u E y E = −u T y where the minus is due to the fact that the power coming in in a system (+ convention) is going out the other system (− convention) whereas the overall balance is zero. Since the robot extracts energy from the tank to implement the desired command u d it is well possible that the total energy can become negative, implying the overall system is no more passive and so possibly unstable. To avoid such situation, we need to introduce the following power preserving interconnection Σ I , 28 as a function of a modulation factor E ∈ [0, 1] where u d is the desired command and u = E u d is its modulated version. Since it is not possible to estimate in continuous-time the energy required by a command but only its actual power u d (t) T y(t), a common solution is to set a threshold T mod . Namely, when T(x E (t)) ≥ T mod no modulation is needed because the tank has enough energy, while, when the level of the energy is smaller than such a threshold, we compute the scaling factor E ∈ [0, 1] according to and derive the corresponding modulated command u = E u d that guarantees the passivity (i.e., T(x E (t)) > T min ). In the modulation case, the value u is as close as possible to u d according to the actual available energy T(x E (t)) − T min . It is easy to see that no command is applied (u = 0) when T(x E (t)) ≤ T min . The threshold T mod is application-dependent as T min , T max : a common way to set them is to record some time-series for command torques and positions over mild working conditions, and derive the time-series for the energy, and finally tune proper values for the thresholds. Unfortunately, no one-size-fits-all recipe can be provided.
It is worth mentioning that the system connection of SPHSs is a SPHS as well. 29 Remark 2. The energy stored in the tank is the mean energy dissipated by the system, implying that T(x E (t)) is deterministic. The negative sign in front of the D S (x) in Equation 6 can be interpreted as a "bias" due to the noise affecting the system. The energy dissipated by the deterministic subsystem through the power > 0 has at least to compensate for the non-passive behavior of the noise represented by the power D S (x) ∶= 1 2 Tr ] . The overall mean energy power T(x E (t)) can be seen as a passivity margin in the sense that the larger is T(x E (t)), the higher is the energy generated from non-passive operations that the system can absorb, still remaining passive. 3,30 The overall model of the SPHS Σ m endowed with the energy tank Σ T is given by Σ m Remark 3. Up to now, the analysis of passivity and convergence has been carried out using the Itô integration. Similar results can be obtained for the Stratonovich integral where W(t) denotes a standard Brownian motion adapted to a reference filtered probability space ( Ω,  , ( t ) t≥0 , P ) and •dW(t) denotes the integration in the sense of Stratonovich. 18,31 Since any Stochastic Differential Equation (SDE) in the Stratonovich sense can be uniquely associated to an SDE in the Itô sense, by introducing the compensation term (⋆) 18 with̃(X(t)) being the n−dimensional vector whose i th component̃i is given bỹ where i,j denoted the entry in position (i, j) of the n × d matrix , see, Reference 18. Thus, the SPHS (11) reads in Itô form as We stress that the choice between the two types of integration is mostly motivated by applications. The geometric nature of PHS would suggest the choice of Stratonovich integral. 29 Nonetheless, since passivity is merely a probabilistic property requiring the use of the infinitesimal generator associated to the driving process, the easiest way to characterize passivity for a system driven by a Stratonovich type of noise, is to reformulate the stochastic system in terms of Itô, hence exploiting the associated martingale property.

Energy tank for stochastic port-Hamiltonian systems with additive noise
In previous subsection, we have considered a PHS with a multiplicative linear noise of Brownian type, where the noise is intended in the sense of Itô. We have mentioned that the choice of a multiplicative linear noise is mainly motivated by stability reasons. In fact, in the deterministic case passivity implies the stability of the autonomous system. A similar result holds within the stochastic framework if the noise vanishes at the stationary solution.
In the present section we consider the case of an additive noise, that is we consider Equation 1 with constant diffusion term, that is, (X(t)) ≡ ∈ R n×d , such that T is a positive definite matrix, where W(t) denotes a standard d−dimensional Brownian motion adapted to a reference filtered probability space . From a purely stochastic point of view, an additive noise is simpler to treat. In fact, in such a setting Itô and Stratonovich stochastic integrals coincide and, integrating by parts, they can be defined in a standard way as a classical integral of a bounded variation function (see Appendix A). Nonetheless, from a stability point of view the constant diffusion complicates the situation under analysis, since the noise does not vanish at the stationary point. This does not allow the process to converge to the desired solution. Therefore, also passivity is a delicate point to study and, as to overcome such an issue, the weaker Definition 2 of ultimate stochastic passivity should be used. 20 Broadly speaking, the passivity of the system is not defined on the whole space but only outside a ball centered at a specified state. This definition has several desirable implication regarding the limiting distribution of the system. It turns out that this weak notion of passivity is tailor-made to deal with SDE characterized by additive noise. Moreover, such notion allows to extend previous results to consider also the case of a SPHS with non vanishing noise.
In the following, we will consider the Hamiltonian H to be a quadratic function of the state, that is, with P a positive definite symmetric constant matrix of suitable dimensions. The following proposition holds. (14) is strictly ultimately stochastic passive.

Proposition 2. The SPHS Σ a in
Proof. The infinitesimal generator of the Itô process SPHS (14)  is From the quadratic form of the Hamiltonian function it follows that Hence, there exists a positive constants C > 0 and > 0 such that, for ||x|| ≥ C it holds and, so for all ||x|| ≥ C, we have that there exists C > 0 so that which is the Definition 2 of ultimately stochastic passivity. ▪ The notion of ultimately stochastic passivity can be thought as follows: for X converging to X * = 0 we have that the noise keeps injecting energy preventing the system from stabilizing at X * = 0. Nonetheless, the system cannot exhibits non stationary behavior since, as soon as the process exits from a ball of radius C, the system becomes passive and the energy injected by the noise into the system is strictly less than the one dissipated, implying that the system recovers stability. Therefore, for large times, the system will converge towards a stationary law.
Remark 4. Proposition 2 highlights why a weaker notion of passivity is necessary. Since P(x) is positive-definite, it is always possible to find x small enough such that inequality (15) is violated.
Since the SPHS (14) is ultimately stochastic passive, we can endow the system Σ a with a (virtual) energy tank as we did for the SPHS with linear multiplicative noise Σ m . We denote the power dissipated  + (x) > 0 by the stochastic system (14), if positive, by Tr where (f ) + ∶= max{f , 0} is the positive part of a given function.
Remark 5. In Equation (16) the positive part of the dissipated energy is stored into the tank. This emphasize a key difference of the proposed weaker notion of passivity. In fact, in the linear multiplicative case (6) we required the validity of an assumption ensuring positivity of the dissipated energy . On the contrary, exploiting a weaker notion of passivity, no positivity assumption has to be made. Intuitively, the new weaker notion of passivity allows the system to have a non-passive behavior without employing any energy stored into the tank to recover passivity. These passive actions are allowed since they happens in a suitable domain, where the energy coming from the noise is more relevant than the true energy of the system. Nonetheless, the general weaker passivity property ensures that the energy due to the noise cannot concur to the overall instability of the system, so that above non-passivity due to the noise is limited to a certain region.
The SPHS endowed with the energy tank is given by where the same considerations reported in Subsection 2.1 hold here.
Remark 6. The notion of strictly ultimately stochastic passivity is fundamental in proving the existence and uniqueness of an invariant measure. In general to ensure stability of a stochastic system it is not enough to require only ultimately stochastic passivity, namely that there exists C > 0, so that for all ||x|| ≥ C, it holds As a counterexample consider the system and by choosing as Lyapunov candidate it can be seen that there exists C > 0 so that, for ||x|| ≥ C it holds H(x) < 0.
Nonetheless, system (18) diverges, hence becoming unstable. To avoid such a phenomenon, thus ensuring system stability, the correct requirement is that there exist C > 0 and > 0, such that for ||x|| ≥ C.

THE STOCHASTIC TELEOPERATION SYSTEM
A bilateral teleoperation system consists of a master robot and a slave robot connected by a communication channel. We consider two n-degree of freedom (n-DOF) and gravity-compensated manipulators where = m, s stands for master/slave, q = (q 1 , … , q n ) ∈ R n is the set of generalized coordinates, C is the Coriolis and centrifugal term, F is the generalized force due to the command, and F ext, denotes the force due to interaction with the external world. From the slave point of view, such force represents the interaction with the environment, hence we indicate it by F s = F e , whereas, from the master point of view, it represents the interaction with the human, then we consequently we indicate it by F m = F h . R m and R s are semi-positive definite matrices modeling viscous frictions present and any additional damping injection. We also included the stochastic noise:̇denotes the formal time derivative of a white noise, with a suitable function representing the noise intensity. The term̃is defined in Equation (12) and it is the correction term needed to consider the noise in Itô sense. In fact, according to the Wong-Zakai approximation theorem, 32 without the suitable correction terms̃Equation (19) should be considered to be in Stratonovich form. Since, we aim at consider the noise in the sense of Itô integration, we added a compensation term. Moreover, we also assume that gravity effects are locally compensated.
The following properties hold for the deterministic subsystem: 33 (i) the matrices M (q ), = s, m are lower and upper bounded, that is, where m (M (q )), resp. M (M (q )), denotes the minimum, resp. maximum, eigenvalue of the matrix M (q ); (ii) the Coriolis matrices C (q ,q ) are given by where Γ j kh (q ) are the Christoffel symbols of the first kind with the symmetric property that Γ j kh (q ) = Γ j hk (q ).
We thus have thatṀ are skew-symmetric matrices, that is,Ṁ (iv) the dissipation matrices D are symmetric and semi-positive definite.
By introducing the generalized momentum p = M (q )q , p = (p 1 , … , p n ) ∈ R n , = m, s, the system (19) can be rewritten using the SPHS formalism 8,34 as where R (q , p ) is the dissipation matrix and is the kinetic energy of the system. Since the robot is gravity-compensated we do not have to consider the potential energy in the definition of the Hamiltonian H (p ). It is worth remarking that, since in Equation 19 we have added a compensation term, the noise in Equation 20 has to be intended in the senso of Itô.
In this article we consider the following two types of noise: the linear multiplicative noise and the additive noise.

Hypothesis 1.
Either assumption (v) ′ or (v) ′′ holds for the stochastic system (20): The assumption (v) ′ refers to the multiplicative noise setting: in such a case, the teleoperated system in free motion (i.e., with F h = 0, F e = 0) converges to the stable solution p = 0 since the noise vanishes for zero velocity. The second assumption (v) ′ refers to the additive noise setting. As to study this case the ultimately passivity must be considered, which leads to a different notion of convergence, related both to the ergodicity and to the invariant measure of the system.
Before stating the main results, we will briefly introduce in Sections 3.1 some key aspects of stochastic delayed systems, while in Section 3.2 we will recall the architecture we will use to control the bilateral teleoperation system.

Stochastic delayed systems
The main problem arising when one aims at obtaining analytical results for stochastic delayed system is that classical methods cannot be directly applied: In particular, the Itô formula cannot be used. To overcome this problem we lift the stochastic process of interest to an infinite dimensional space. Such procedure allows to exploit tools derived within the theory of stochastic analysis in infinite dimension. 35 A standard method to obtain the aforementioned lift consists in introducing the so called segment process 35 and taking values in the space of continuous functions  ∶=  ([− , 0]). Adding an arbitrary initial condition ∈  to the original system, the existence and uniqueness of a solution to the initial value Cauchy problem in  can be thus proved.
A formal definition of the Itô formula for infinite dimensional processes has been a long-standing problem, mainly due to severe difficulties in finding a proper definition for the quadratic covariation process. During last years some results appeared proving formally Itô formula for infinite dimensional stochastic process. 35,36 A particular attention has been given to Itô formula for stochastic delay equations, that has been proven in several settings, see, for example, References 35,37. To avoid the general mathematical machinery needed to define an infinite dimensional Itô formula, we can exploit the particular form for the considered delay in (21). Such delay is usually mathematically called discrete delay. More precisely, a discrete time-delay can be defined as ∫ 0 q(t − ) (d ) with represented by a sum of Dirac measures, namely ∶= ∑ x i ∈(0, ] x i , where x i is the Dirac measure centered at x i . Let us recall the following trick typically exploited to formally obtain the Itô formula, for such type of delay. Let , it can be seen that for t ∈ [0, ), the effect of the delay acts as a parameter since the delayed value of the system X(t − (t)) is known using the initial condition . The main advantage of this reasoning is that standard results from (finite dimensional) stochastic analysis can be used, with particular reference to the Itô formula. On the next time interval [ , 2 ) an analogous reasoning holds, since the delayed-valued of the system has been evaluated at the previous step; in particular the time delay t − (t) ∈ [0, ) and the corresponding value of the process X(t − (t)) is known from previous step. By iteratively constructing the solution, we can notice that, for any interval of the form [n , (n + 1) ), n ∈ N, the solution can be derived using known values taken within previous subinterval [(n − 1) , n ), and again the standard finite dimensional Itô formula holds.

F I G U R E 1
The two-layer architecture: passivity layer in blue and transparency layer in green. The operator and the environment are modeled as impedance systems mapping velocity into (generalized) force,q m  → F h andq s  → F e ; The master and slave robots have an admittance causality mapping force into velocity/position, F h − F m  → (q m ,q m ) and F e − F s  → (q s ,q s ); The controllers in the transparency layer at the master and slave side are PD with parameters K and B as in Equation (22); m and s are the communication delays from master-to-slave and slave-to-master, respectively; T m and T s are the energy tanks within the passivity layer that modulates the commands if the energy is not enough to guarantee passivity-preserving commands, that is, F m = E m F dm , F s = E s F ds

Two-layer teleoperation architecture
In the present work we will consider the bilateral teleoperation architecture where the coupling between the master and slave robots is implemented via a two-layer approach. 3,10 In the first layer, called transparency layer, the desired coupling forces, namely F dm and F ds , are computed according to any kind of controllers. Such control forces go through the passivity layer implementing the commands using only the energy available within the tank: If not enough energy is stored in the tank, the passivity layer modulates the force command to be implemented. Moreover, if the energy stored in the tank is below a certain threshold, no energy extraction from the energy tank is allowed: no coupling force is thus implemented by the passivity layer as explained in (8) and (9). Figure 1 shows the block diagram of the two-layer architecture where we implemented in the transparency layer a PD coupling, that is, spring-damper coupling, 38 both at the master and at the slave side ) .
Matrices K and B are positive definite, whereas m (t) > 0 and s (t) > 0 are the master-to-slave and slave-to-master time-varying delays. 8 It follows that the master and slave robots exchange position and velocity information. We will show that energy tanks T m and T s ensure passivity even if a delay is present also in the stochastic settings. Due to communication delays, in a bilateral teleoperation system we need two tanks: one on the master side, T m , and another on the slave side, T s . Since passivity is a property of the whole system (master side + communication channel + slave side), it is advisable to allow the exchange of energy between the two tanks to avoid the situation where one of the tanks is close to the minimum threshold T min , while the other has more energy than necessary. Therefore, the tank model described in (7) has to be enhanced to take into account the energy flow between the master and the slave side as shown in Figure 1. Each tank will then be modeled as where if = 1 the tank stores the dissipated energy  (x) and the incoming power P in (t), whilst it releases the outgoing power flow P out (t) if there has enough energy (i.e., T(x E (t)) > T ava ), if = 0, otherwise. The energy flows are modeled as For example, if the energy stored in the master tank T m goes below a fixed threshold T m req , then an energy request E m req is sent to the slave tank. If the energy stored in the slave tank T s is greater than a fixed threshold T s ava , then the needed energy is sent back to the master tank. An analogous reasoning is valid for an energy request sent by the slave tank.
The application-dependent constant P is the rate of energy flowing from one tank to the other, whereas the constants in (24) satisfy the following inequalities where the meaning of lower limit T min and the upper limit T min has been explained in Section 2.1. We refer the interested reader to the literature, see, for example, References 5 and 10, for a more detailed explanation of this mechanism.
Remark 7. It is worth recalling that passivity is a sufficient condition for stability. The tank is the mathematical tool that guarantees the passivity even though non passive elements, as the communication delay, could destabilize the overall system. A non passive action that could bring to instability should drive the energy within the tank to go below zero (the system is no longer passive but it is active and so unstable). The tank modulation mechanism explained in (8) and (9) guarantees that commands that could bring the energy below the minimum threshold are modulate and, ultimately, zeroed. In this way the system is passive by construction (energy strictly larger than zero). The price to pay for the passivity (and so stability) is performance degradation. Anytime the passivity layer modulates the commands F dm and F ds computed by the transparency layer the performance perceived by the operator worsen. However this is the mechanism exploited by the tanks to preserve stability. About the choice of the upper and lower bounds, it is impossible to provide a general recipe since they are strongly application-dependent. The thresholds used in a teleoperated surgical robots are order of magnitude smaller that threshold used in a teleoperator moving heavy tools. Moreover force/velocity scaling comes often into play to improve the usability of a teleoperation system as shown in for example, Reference 39. The proper choice of the thresholds involve a tuning similar to the one needed to choose the controller gains K and B.
In the following two sections, we will tailor this bilateral teleoperation architecture for the cases of SPHSs with multiplicative and additive noise.

THE STOCHASTIC TELEOPERATION SYSTEM WITH MULTIPLICATIVE NOISE
In this section, we consider the master-slave SPHSs with multiplicative noise endowed with two energy tanks of the form where  (p ) is given in Equation (6). Explicit calculations yield It therefore follows that where Λ ∶= T is related to the variance of the multiplicative noise.
Remark 8. In Equation (27) a slight abuse of notation has been used. In fact, the dissipated energy stored into the tank should be formally equal to the value of the process p conditioned to the filtration  t , that is, , t ≥ s. (27) is motivated to avoid unnecessary heavy notation.

The choice made in Equation
The coupling forces are implemented using the energy stored in the tanks T m , T s . The power port of the tank (u E , y E ) is interconnected with the power port of the robot (F ,q ) through the power preserving interconnection in (8). In particular, if F T dq < 0, then, as to implement the input F d , some energy must be extracted from the tank; on the contrary if F T dq > 0 then the input is dissipative, therefore it increases the available energy of the overall system.
The following two results show the stochastic passivity of the teleoperation system and the asymptotic stability when the system has no exogenous inputs.

Proposition 3 (Stochastic passivity).
The master-slave stochastic bilateral teleoperation system described by two SPHS with energy tanks (26) and spring-damper coupling (22) is stochastically passive at the port Proof. To prove the Proposition, we exploit the technique introduced at the beginning of the previous Section, hence by discretizing the time-interval as to apply the Itô formula. Consider the total energy of the system (26) and denote ∶= min{ m , s }. For t ∈ [0, ), let be the initial condition for the delayed system, that is, The standard Itô formula for t ∈ [0, ] yields to Since p , = m, s, and using the assumption R M −1 − 1 2 Λ ≻ 0, it follows that  > 0. Therefore, for t ∈ [0, ], we have Iterating for any t ∈ [n , (n + 1) ), n ∈ N, as explained above, we obtain the claim for any t ∈ [0, ∞). ▪ In the following we will denote ) .
We have the following convergence result.

Proposition 4 (Asymptotic stability).
The master-slave stochastic bilateral teleoperation system described by two autonomous SPHS (i.e., F h = 0 and F e = 0) with energy tanks (26) and spring-damper coupling (22) is asymptotically stable in probability. In particular, the system reaches a full stop and the position tracking is achieved, that is, Proof. Using the passivity property proved in Proposition 3 in free motion, that is, F h = 0 and F e = 0, we have that H(X(t)) ≤ 0, implying that the process H(X(t)) is a super-martingale. Using Dynkin-formula it follows that H(X(s))ds < 0, so that, from Chebyshev's inequality, we obtain from the super-martingale property that for > 0 it holds with c a constant depending on > 0. Therefore X converges in probability. 40(theorem 2.2) .
Letting t → ∞ we obtain that (p m , p s ) converges in probability toward (p * m , p * s ) = (0, 0). Therefore, since the momentum p converges to 0 as t → ∞, the system reaches a full stop.
From Reference 41 (lemma 5.3) together with the fact that super-martingale admits almost surely a long-time limit, it follows lim t→∞ |p (t)| = 0.
We therefore have that Using thus the system dynamics (26) with the spring-damper coupling (22), the asymptotic convergence in probability to 0 of p also implies the position tracking (29). ▪ Remark 9. An argument analogous to the one in Remark 3 shows that Proposition 3 holds true in the case of Stratonovich multiplicative linear noise if the following condition is valid In fact, taking into account the Stratonovich correction term as in Remark 3, it follows that implying that if (30) holds than the system is passive.

THE STOCHASTIC TELEOPERATION SYSTEM WITH ADDITIVE NOISE
We now consider the master-slave stochastic PHS, for = m, s, with additive noise endowed with two energy tanks where  + is defined in Equation (16) and it explicitly reads in the present situation as Considering the same coupling between tanks and systems, and the spring-damper coupling as controllers, we can prove the results that correspond to Propositions 3 and 4.
Proposition 5 (Ultimately stochastic passivity). The master-slave stochastic bilateral teleoperation system described by two SPHS with energy tanks (31) and spring-damper coupling (22) is strictly ultimately stochastic passive at the port ) .
Proof. Proceeding as in the proof of Proposition 3, we iteratively apply the Itô formula to any t ∈ [n , (n + 1) ], with ∶= min{ m , s }. Thus, for t ∈ [n , (n + 1) ], it holds Using Equation (16) we have that for a suitable constant C > 0, it holds for ||p || ≥ C, for = m, s and for > 0. Therefore, we have that there exists C > 0 such that for ||p || ≥ C, obtaining the claim. ▪ In order to prove the convergence of the SPHS toward an invariant measure, the Feller property and the transition Markov semigroup for the SPHS have to be studied. Since the control commands (22) are delayed due to the communication channel, to show the Markov property of the process, it must be lifted in a suitable infinite dimensional space. We will first recall basic definitions and results about ergodicity and Markov property for stochastic delay equations, referring the reader to References 21,22,42 for further details.
Introducing the segment process as done in Section 2, it can be shown, 22(proposition 3.3) that the SPHS (31) is a Markov process in the space  ∶= ([− , 0]). Recall that the space  is a separable Banach space when equipped with the norm || ⋅ || ∞ defined as As in Section 4, we will denote for short ) .
Further, the notation X t ∈  denotes the segment process, defined as It can be therefore shown that the segment process (X t ) t≥0 is a Markov process on , in the sense that for all t ≥ s and Borel set B ∈ (), where () denotes all the Borel set over ; see Reference 22 (proposition 3.3). Since the Markov property holds in the space , the Markov transition function can be introduced, namely where X t is the infinite dimensional lifted segment-process X with initial value ∈ . Following Reference 42 (chap.9), we will say that the transition semigroup p is called Feller semigroup if the Markov semigroup defined as is bounded and continuous for any H ∈ C b (). If P t H( ) is continuous and bounded for any t > 0 and for any H ∈ C b (), then it is called strongly Feller semigroup. If, for t > 0, all Markov semigroup P t are equivalent, then P t is called t−regular.
The Markov semigroup and the transition function are connected by the following (t, , B), B ∈ ()).
For the definitions of strongly Feller Markov semigroup and regular Markov semigroup, that will be used later, we refer to the literature. 21,42 For the sake of readability, before stating and proving the main results, we highlight the main steps to prove existence and uniqueness of an invariant measure. Regarding uniqueness, we will prove that associated transition semigroup is strongly Feller and regular. In order to do that, we exploit the fact the Feller property is invariant under change of measure, 26(theorem 2.1) , so that we can reduce the original delay equation to a standard stochastic equation affected by no delay. A similar argument has been used in References 22,43. Since, we reduced the system to a finite dimensional stochastic process, which is Feller by standard arguments, we have only to prove that the Feller is indeed inherited by the segment process. This at last implies the uniqueness of any stationary solution. For what concern existence of an invariant measure, we make use of the ultimately stochastic passivity property of the system, Proposition 2, to show that there exists an invariant measure. In particular, the existence follows from the next result, 44(proposition 3.1) , that we report in order to make the treatment as much as self-contained as possible. Proposition 6. ( 44 proposition 3.1) Let (X, || ⋅ || X ) be a complete separable metric space and let (P t ) t≥0 be the semigroup of Markov operators corresponding to a Markov process which satisfies the Feller property. Assume that there exists a compact set B and a point x ∈ X such that Then the Markov process has a stationary distribution.
Next Proposition proves a convergence result for the stochastic teleoperated system with additive noise.
Proposition 7 (Asymptotic stability of the transition measure). For the master-slave stochastic bilateral teleoperation system described by two autonomous SPHS (i.e., F h = 0 and F e = 0) with energy tanks (31) and spring-damper coupling (22), there exists a unique invariant measure . Furthermore the following holds true as t → ∞.
To prove Proposition 7, we need a complementary lemma.
Proof. We start proving the strong Feller property. It is worth remarking that the strong Feller property is invariant under the change of measure, 26(theorem 2.1) . Applying Girsanov Theorem 45 to q , we can introduceW , a Brownian motion under the probability measureP equivalent to the original probability measure P defined as The overall dynamics for p (i.e., robot, controller, tank, and communication channel) can be rewritten as dp m (t) = f m (q m , p m )dt + m dW m (t), dp s (t) = f s (q s , p s )dt + s dW s (t), with Notice that, after application of Girsanov theorem, no delay affects the system (33). Then the transition probabilityp(t, x, B) of the new system under the measureP is equivalent to the original transition probabilityp(t, x, B) under P. 46 Let us also underline that, as above, in order to apply standard Girsanov theorem we must consider the sequence of interval [n , (n + 1) ).
Let us underline that a delay term appears in Equation 34 through P in . This implies that the Equation 34 continues to be a delayed equation. Nonetheless, no delay affects q and p . Therefore, standard arguments belonging to the theory of finite dimensional Markov processes can be exploited.
Let H be a bounded functional on . Then, we have to prove that the map is continuous for ∈ . Notice that, being (q , p ), = m, s, within non-delayed Equation 34, they only depend on (0). Moreover, it is well-known that hypothesis 1 implies that the couple (q , p ) in (34) generates a strongly Feller semigroup, see Reference 21 (theorem 7.1.1). Therefore, further using the fact that x E continuously depends on the couple (q , p ), = m, s, and exploiting both the Itô formula and the Gronwall lemma, we have that, for 1 , 2 ∈ , it holds then, by Reference 23 (theorem 1.7) and for t > , it holds allowing us to conclude that the segment process X is strongly Feller continuous for t > . So that, by Reference 22 (proposition 5.2), we infer the continuity in  of the map (36) and the Feller property holds, for t > .
Following Reference 47 (corollary VIII.2.3), see also Reference 22 (corollary 5.3), the process X in Equation (35) is irreducible at time t − , which together with the strong Feller property implies 21(proposition 4.1.1) the regularity for t > 2 .▪ We are now ready to prove Proposition 7.
Proof of Proposition 7. (Existence). A successive application of Itô formula for any interval [n , (n + 1) ) yields following Proposition 5 that in the case of free-motion, that is (F h , F e ) ≡ 0, there exists > 0 such that for ||x|| ≥ C, C > 0. Then, Itô formula implies Using Equations (37) and (38) being C m the maximum value attained by  over the set |x| ≤ C.
Using the fact that (x) ≥ 0, it follows from Equation (38) E(X(s))ds, so that using estimate (39) we obtain The existence of the invariant measure thus follows using Equation (40)

NUMERICAL IMPLEMENTATION
The current section aims at illustrating the results derived along previous Sections both in the case of a linear multiplicative noise and in the additive noise case. The implemented numerical examples will show how different types of noise imply different notions of stochastic stabilization. We focus only on the free motion scenario. Namely, we analyse the case when the slave robot does not interact with the environment, F e = 0, as to show how the trajectory error behaves for the two different noises. Moreover, we assume that F h is an exogenous signal independent of the velocity of the master robot.

Linear multiplicative noise
Consider a 2−dimensional master-slave system as given in Equation 26. For the sake of simplicity we will consider two second-order linear time-invariant systems for the master and slave robots with coefficients given by   . Shaded region is when the external force F h affects the system. From Figure 2 it emerges how, under the proposed control scheme equipped with an energy tank, position tracking is achieved in free-motion. Figure 2 (bottom panel) shows the level of the energy tank in the multiplicative noise case.

Additive noise
We now consider the stochastic system affected by additive noise as given in Equation 31, with the same parameters as in Section 6.1. Figure 3 reports . Shaded region shows when the external force F h affects the teleoperation system. Differently from the multiplicative noise case, the system is not asymptotically stable due to a non vanishing noise, Figure 2. Nonetheless, it can be seen how the overall system is stable in the sense that the system reaches a limiting invariant distribution; as a consequence error tracking does not converge to the equilibrium point (q m = q s ) but kept oscillating around it. Also, it is worth highlighting how the tracking error strongly depends on the noise: a smaller noise along the y−axis implies a smaller oscillation in the position error, whereas a larger noise along the x−axis implies a larger oscillation. This can be easily seen in the steady-state regime (t > 55 s). Figure 3 (bottom panel) shows the level of the energy tank in the additive noise case.

CONCLUSION
In this article we have shown how a stochastic master-slave bilateral teleoperation system can be endowed with energy tanks as to ensuring the overall passivity of the system despite unknown environment and unknown and time-varying communication delays. The proposed setting presents several novelty when compared to the deterministic case. In particular, it allows the noise to stochastically inject energy into the system, forcing tanks to take it into account. Accordingly, we have shown how the energy tanks can be efficiently used to force passivity when the master and slave robots are interconnected via a PD-like control. First addressing the case of a stochastic system with linear multiplicative noise, then considering the additive case, we have shown how to endow the system with an energy tank in such a way that the stored energy can be exploited to passify the overall teleoperated system. It is worth underlying that the additive case constitutes the most delicate one when the final goal is to prove passivity or stability of the system. Therefore, a non-vanishing noise has to be imposed to prevent the system from being neither passive nor stable. Moreover, as to also tackle a non-vanishing noise setting, an ad hoc weak notion of passivity has been introduced, allowing us to derive an associated weak notion of convergence for the corresponding stochastic system. We remark that the introduction of the above mentioned weaker notion of passivity turn out to be strictly necessary when an additive noise perturbs the system since, in this scenario, the system fails to be passive, or stable, in any classical stochastic sense.

APPENDIX A. ON THE CHOICE OF STOCHASTIC INTEGRATION
The present section is devoted to a brief review of stochastic integration with the goal of highlighting the main differences of the possible notions of stochastic integrations. We refer to References 45 and 18 for a detailed introduction on stochastic integration.
Equation (A6) is often written using the short-hand notation dX(t) = (t, X(t))dt + (t, X(t))•dW(t), for a Stratonovich stochastic differential equation (SDE), or dX(t) = (t, X(t))dt + (t, X(t))dW(t), for a Itô SDE. The Itô integral ∫ t 0 (t, X)dW(s) does not look into the future", as in the definition of Itô integral the left point of the partition is taken, that is, t * j = t j ; this implies the remarkable property that the Itô integral is a martingale whereas the Stratonovich integral is not. In particular this fact allows to interpret the coefficient in Equation A8 as the mean change in X, usually referred to as drift, whereas is the volatility representing the randomness in the system; similar interpretation does not hold for the corresponding version in Stratonovich sense. The martingale property is the main reason why in biology or in finance the Itô integral is the most common choice.
However, the Itô integral does not satisfy the chain rule of calculus, so that classical rule of differential calculus cannot be applied. In particular, in proving the chain rule, the unbounded variation property of the Brownian motion implies that the second order term in the Taylor expansion does not vanish. More precisely, given a general function f , then it holds df (X(t)) = f ′ (X(t))dX(t) + 1 2 f ′′ (X(t)) 2 (X(t))dt, being X the solution to the Itô SDE (A8). On the contrary, if we have assumed that X is the solution to the Stratonovich SDE (A6) the classical chain rule of calculus would have applied, that is, df (X(t)) = f ′ (X(t))•dX(t), see, for example, Reference 18. Further, mostly in physical applications, differential equation (DE) has to be properly defined on a general manifold. In particular requiring that the coefficients of the DE belongs to the tangent space ensures that the solution of the DE stays on the manifold. Unfortunately, also this does not immediately generalize to the stochastic case. In particular, given  a general manifold, it can be shown that requiring that and ∈ T, being T the tangent space to , then X lives on the manifold  if the integration is intended in the Stratonovich sense, that is X is the solution to the SDE (A6). The same does not hold if the integration is intended in the Itô sense. This last two properties are at the core why in physics, and in particular in geometrical mechanics, the choice of Stratonovich integration is usually made.
In a very broad generalization, the Stratonovich integral enjoys good geometrical properties whereas the Itô integral has good probabilistic properties; therefore the physics and engineering community usually consider SDE driven by a Stratonovich noise whereas researchers in biology and mathematical finance consider the Itô integral.
At last it has to be mentioned that the two integral interpretations are related. The solution X to the SDE in Stratonovich sense in Equation (A6) can be rewritten in Itô form as (t, X) + ′ (t, X) (t, X) ) ds + ∫ t 0 (t, X)dW(s), see, for example, Reference 18. This last property suggests that typically, the definition of the integral that is most suited to the application is decided at the beginning, and then one can pass from one notion of stochastic integral to the other, recovering the definition that allows to obtain more easily the needed results. For instance, as done above, being the very definition of PHS of geometric nature, the most natural notion of integral to consider is in Stratonovich sense, as done, for example, in Reference 29. Nonetheless, when aiming at defining passivity properties of the system, the most simple notion of integral to consider is the Itô version, since the martingale property can be exploited to easily carry out all the computations; therefore the Stratonovich SDE is converted into the corresponding Itô SDE.
However, let us stress that, although it is possible to convert a Stratonovich SDE into a Itô one, the choice of the integration has to be made at the beginning. In fact different notions lead to different results; for instance, it is shown in