Boundary feedback control of an anti-stable wave equation
BBOUNDARY FEEDBACK CONTROL OF AN ANTI-STABLE WAVEEQUATION
PIERRE APKARIAN AND DOMINIKUS NOLL Abstract.
We discuss boundary control of a wave equation with a non-linear anti-damping boundary condition. We design structured finite-dimensional H ∞ -output feed-back controllers which stabilize the infinite dimensional system exponentially in closedloop. The method is applied to control torsional vibrations in drilling systems with thegoal to avoid slip-stick. Key words:
Wave equation · boundary feedback control · anti-damping boundary · infinite-dimensional H ∞ -control · slip-stick · torsional vibrations · large magnitude sectornon-linearity Introduction
We discuss H ∞ -boundary feedback control of a wave equation with instability causedby boundary anti-damping. This is applied to the control of vibrations in drilling devices.The system we consider is of the form G nl : x tt ( ξ, t ) = x ξξ ( ξ, t ) − λx t ( ξ, t ) 0 < ξ < , t ≥ x ξ (1 , t ) = − x t (1 , t ) + u ( t ) αx tt (0 , t ) = x ξ (0 , t ) + qx t (0 , t ) + ψ ( x t (0 , t )) (1)where ( x, x t ) is the state, u ( t ) the control, and where the measured outputs are y ( t ) = x t (0 , t ) , y ( t ) = x t (1 , t ) . (2)The non-linearity ψ satisfies ψ (0) = 0 , ψ (cid:48) (0) = 0 , and the steady state is ( x, x t , u ) =(0 , , . The linearized system G is obtained from (1) by dropping the term ψ ( x t (0 , t )) .The parameters satisfy λ ≥ , α ≥ , while q is signed. System (1) was first discussed in[46, 18, 47] in the context of oil-well drilling. The author of [18] proves open-loop stabilityof (1) for the case q < using Lyapunov’s direct method. Since applications typicallylead to the opposite case q > , where instability occurs, various control strategies havebeen proposed for that setting.Lyapunov’s direct method is used in [39, 40, 12, 13]. This leads to infinite-dimensionalcontrollers which are not implementable, or to observer-based controllers, which due totheir lack of robustness are out of favor since the late 1970s.Delay system techniques are used in [24, 38, 39, 41, 21, 4], but require λ = 0 , whichleads to an oversimplified model. Input shaping is used in [32], but as presented, alsorequires the un-damped model λ = 0 .Backstepping control is used in [43, 33, 17, 34, 35, 23], but with the exception of [17],where λ = α = 0 , leads to infinite dimensional or state feedback controllers, which are notimplementable. Infinite dimensional controllers can also be obtained with the method in[11]. Other ideas to avoid slip-stick include the design of feedforward startup trajectories[1], or manipulation of the weight on bit in [38, 39]. Model (1), (2) has also been used tocontrol axial vibrations, see [39, 14], and for robotic drilling [14]. ONERA, Department of System Dynamics, Toulouse, France. Institut de Mathématiques, Université de Toulouse, France. a r X i v : . [ m a t h . O C ] F e b PIERRE APKARIAN AND DOMINIKUS NOLL What these approaches have in common is that they are guided by the method ofproof of infinite-dimensional stability. This leads to impractical control laws. In contrast,approaches guided by practical considerations have also been applied to oil-well drilling[42, 15], but those use finite-dimensional approximate models. This makes it desirableto bridge between both approaches by designing practical controllers using the infinite-dimensional model (1). In the present work we design H ∞ -controllers with the followingrequirements:(a) The controller is output feedback and of a simple, implementable structure, like areduced-order controller or a PID.(b) The controller stabilizes the infinite-dimensional system (and not just a finite-dimensional approximation of it).(c) H ∞ -optimality of the controller is certified in closed loop with the infinite-dimensionalsystem (and not just with a finite-dimensional approximation).(d) Due to the achieved infinite-dimensional H ∞ -performance, slip-stick is avoided, orat least mitigated.These requirements are achieved by going through the steps of the following general H ∞ -control scheme, which we proposed for boundary and distributed control of PDEs in[9, 5, 4], where it has already been applied successfully to a variety of applications. Algorithm 1.
Infinite-dimensional H ∞ -design (cid:46) Step 1 ( Steady-state ) . Compute steady state of non-linear system G nl and ob-tain linearization G . Compute transfer function G ( s ) and determine number n p ofunstable poles of G . (cid:46) Step 2 ( Stabilize ) . Fix practical controller structure K ( x ) , and compute initialstabilizing controller K ( x ) for G . Use Nyquist test to certify stability of linearinfinite-dimensional closed loop. (cid:46) Step 3 ( Performance ) . Determine plant P with H ∞ -performance and robustnessspecifications, addressing in particular the non-linearity. (cid:46) Step 4 ( Optimize ) . Solve discretized infinite-dimensional multi-objective H ∞ -optimization program using a non-smooth trust region or bundle method [5, 9]. (cid:46) Step 5 ( Certificate ) . Certify final result in infinite-dimensional system withinpre-specified tolerance level as in [5, 9].While some of the elements of algorithm 1 are standard, others need to be adapted tothe current case and to be explained in detail. In section 2, the mechanical model forcontrol G nl will be derived. Its linearization G , transfer function, and open-loop propertieswill be discussed in sections 3 and 4. Locally exponentially stabilizing controllers will besynthesized in section 5, and H ∞ -synthesis for the full, non-linear model in section 6 willcomplete the procedure. Numerical results are regrouped in section 7.2. Model of drilling system
We derive model (1) from the setup of an oil-well drilling system, shown schematicallyin Fig. 1. The state of the system is described by the angular position θ ( ξ, t ) and angularspeed θ t ( ξ, t ) of the drillstring, where position ξ = 0 refers to the rotary table (top), while ξ = L represents the drill bit (bottom hole assembly), with L the length of the string. ONTROL OF AN ANTI-STABLE WAVE 3
The dynamic equation and boundary conditions are
GJ θ ξξ ( ξ, t ) = Iθ tt ( ξ, t ) + βθ t ( ξ, t ) 0 < ξ < L, t ≥ I B θ tt ( L, t ) = − GJ θ ξ ( L, t ) − φ ( θ t ( L, t )) θ t (0 , t ) = GJc a θ ξ (0 , t ) + Ω( t ) (3)where G is the angular shear modulus, J is the geometrical moment of intertia, I isthe inertia of the string, I B is the lumped inertia of the bottom hole assembly, c a isrelated to the local torsion of the drillstring, Ω( t ) is the time-dependent rotational velocitycoming from the rotary table at the top, used to drive and control the system, while theundriven bottom extremity (bit) is subject to a torque φ ( θ t ) representing rock-bit andmud-bit interaction of the drill bit, depending non-linearly on the rotary speed θ t ( L, t ) at the bottom; [18, 45, 46, 47]. The torsional excitations of the drillstring caused by thefrictional force φ ( θ t ( L, t )) lead to twisting of the string, and this effect propagates alongthe structure from bit to top as a wave with damping factor β > . Similarly, alterationsin the rotary speed Ω( t ) at the top are transmitted to bottom by the same damped wave.This implies that a control action at the top will be delayed by one period of the wavebefore taking effect at the bottom. If measurements are taken only at the top, then thedelay before a control action takes effect is even two periods. Introduction
Problem statement
Oilwell drillstrings are mechanisms that play a key role in the petroleum extractionindustry. Failures in drillstrings can be signiÖcant in the total cost of theperforation process. Vertical oilwell drilling system.
S. Mondie & B. Saldivar (CINVESTAV)
Drilling system
October 2012 4 / 39
Drilling Dynamic Dysfunctions All these degrees of freedom in motion allow the development of specific drilling dynamic dysfunctions when the drillstring follows its intention and is rotated, see Figure 2. In the subsequent subchapters three of them will be discussed in more detail. Stick-slip is the one that is primarily focused on throughout this work.
BendingWhirlLateral Vibration Stick SlipTorque ShockTorsional Vibration Bit BounceAxial Vibration Bit Efficiency LossBendingWhirlLateral Vibration BendingWhirlLateral Vibration Stick SlipTorque ShockTorsional Vibration Stick SlipTorque ShockTorsional Vibration Bit BounceAxial Vibration Bit Efficiency LossBit BounceAxial Vibration Bit Efficiency Loss
Figure 2 –
The three main modes of vibration a BHA can be subject to. [4]
Stick-slip is a phenomenon whose occurrence is enabled by the low torsional stiffness of drillstrings. At the surface a drillstring is driven with a constant rotary speed. “
However, the rotary speed at the opposite end of the drillstring, at the bit, oscillates around the surface RPM ” [5] . The RPM oscillations can reach severity levels where the bit comes to a complete stop for a short moment (stick). Due to the continuing surface drive, after the short stick period the bit is forced to catch up the developed downhole to surface revolutions difference. The consequence is a phase of rotational acceleration up to peak velocities of two or three times the surface RPM. As the drillstring is slipping its rotational restriction, this phase is named slip. Sequences of stick and slip phases are known as stick-slip, see Figure 3. Figure 1.
The goal of the active control scenario is to maintain the system at steady state withconstant rotational velocity θ t ( L, t ) = Ω at the drill bit position ξ = L by acting on thedriving rotary force Ω( t ) , and using measurements of the rotary speed at top and bottom.The steady state solution of (3) is easily obtained as θ ( ξ, t ) = Ω t − (cid:18) φ (Ω) + β Ω LGJ (cid:19) ξ + β Ω2 GJ ξ and this corresponds to applying a constant control torque Ω( t ) = Ω at the top, where Ω = Ω + φ (Ω) + β Ω Lc a . Writing the state in the form θ ( ξ, t ) = θ ( ξ, t ) + ϑ ( ξ, t ) for an off-set variable ϑ ( ξ, t ) , andsubtracting the steady state from (3), we obtain the equivalent system GJ ϑ ξξ ( ξ, t ) = Iϑ tt ( ξ, t ) + βϑ t ( ξ, t ) I B ϑ tt ( L, t ) = − GJ ϑ ξ ( L, t ) + φ (Ω) − φ (Ω + ϑ t ( L, t )) GJ ϑ ξ (0 , t ) = c a ( ϑ t (0 , t ) − Ω( t ) + Ω ) (4) PIERRE APKARIAN AND DOMINIKUS NOLL A dimensionless system is now obtained by the change of variables ξ = L (1 − ζ ) τ = 1 L (cid:114) GJI t.
On putting x ( ζ, τ ) = ϑ ( ξ, t ) , this leads to the following equivalent dimensionless form x ζζ ( ζ, τ ) = x ττ ( ζ, τ ) + βL √ GJ I x τ ( ζ, τ ) I B LI x ττ (0 , τ ) = x ζ (0 , τ ) + LGJ (cid:32) φ (Ω) − φ (cid:32) Ω + 1 L (cid:114) GJI x τ (0 , τ ) (cid:33)(cid:33) x ζ (1 , τ ) = − c a √ GJ I x τ (1 , τ ) + c a LGJ (Ω( t ) − Ω ) (5)We re-write the second boundary condition of (5) at ζ = 1 as x ζ (1 , τ ) + x τ (1 , τ ) = (cid:18) − c a √ GJ I (cid:19) x τ (1 , τ ) + c a LGJ (Ω( τ ) − Ω ) . Taking into consideration that the measured outputs of (3) are the angular velocitiesat the top and bottom positions y ( t ) = θ t ( L, t ) , y ( t ) = θ t (0 , t ) , the outputs of thecentered system (5) may be understood as measurements of the offset angular velocities y ( τ ) = x τ (0 , τ ) and y ( τ ) = x τ (1 , τ ) . This allows us to introduce the control u ( τ ) = (cid:18) − c a √ GJ I (cid:19) x τ (1 , τ ) + c a LGJ (Ω( τ ) − Ω ) , which when chosen in feedback form u ( τ ) = K ( y ( τ ) , y ( τ )) leads to the following finalfeedback control law for (3): Ω( t ) = Ω + (cid:20) K ( y ( t ) , y ( t )) + (cid:18) c a √ GJ I − (cid:19) y ( t ) (cid:21) GJc a L , which is linear as soon as u = Ky is a linear controller. With that the second boundarycondition takes indeed the form x ζ (1 , τ ) + x τ (1 , τ ) = u ( τ ) in (1).Switching back for convenience to t for time and ξ ∈ [0 , for the spatial variable, andintroducing the dimension free parameters(6) α = I B LI , λ = βL √ GJ I , q = − φ (cid:48) (Ω) √ GJ I , system (5) turns into the form (1) with the non-linearity given by(7) ψ ( ω ) = LGJ (cid:32) φ (Ω) − φ (cid:32) Ω + 1 L (cid:114) GJI ω (cid:33)(cid:33) − q · ω, with ψ (cid:48) ( ω ) = − φ (cid:48) (cid:16) Ω + L (cid:113) GJI ω (cid:17) √ GJ I − q, ψ (cid:48)(cid:48) ( ω ) = − φ (cid:48)(cid:48) (cid:16) Ω + L (cid:113) GJI ω (cid:17) LI .
From the definition of q it can be readily seen that ψ (0) = 0 and ψ (cid:48) (0) = 0 , which complieswith the requirement in (1). For later use, we introduce an additional parameter p := ψ (cid:48)(cid:48) (0) , which represents the curvature of the non-linearity at the reference position, and givesinformation on its severity. ONTROL OF AN ANTI-STABLE WAVE 5 Ω angular velocity f r i c t i o n f o r ce φ Ω angular velocity f r i c t i o n f o r ce φ Figure 1.
The case (⌦) > (left) leads to a stable open loop [ ? ]. Thepotentially unstable case (right) is when increasing rotary speed reducesfriction. q α Figure 2.
Five regions and five scenarios. The pattern of the gray zone is0. Blue: 0-2-0. Red: 1-2-0. Magenta: 1-0-2-0. Green: 1-0. ψ ( · ) G y K u y W u z u W z Figure 2.
The case φ (cid:48) (Ω) > (left) leads to a stable open loop [18]. Thepotentially unstable case (right) is when increasing rotary speed reducesfriction.Phenomenological models of the frictional force φ ( · ) have been proposed in the litera-ture. For instance [27] considers a model of the form φ ( θ t ) = φ mud ( θ t ) + φ rock ( θ t ) , where the mud friction is assumed of viscous form φ mud ( θ t ) = c b · θ t , while the rock-bitinteraction is the non-linear(8) φ rock ( θ t ) = W ob R b (cid:104) µ cb + ( µ sb − µ cb ) e − γbνf | θ t | (cid:105) sign( θ t ) , φ mud ( θ t ) = c b · θ t . Here W ob is the weight on bit, R b is the radius of the drill, the non-linear term features thestatic and Coulomb friction coefficients µ sb , µ cb ∈ (0 , , while the coefficient γ b ∈ (0 , is the velocity decrease rate accounting for the Stribeck effect. The fact that µ sb > µ cb is the ultimate reason why the slip-stick phenomenon may occur. Namely, for Ω > wehave φ (cid:48) (Ω) = c b − W ob R b ( γ b /ν f ) ( µ sb − µ cb ) e − ( γ b /ν f )Ω , which leads to q = − c b + W ob R b ( γ b /ν f ) ( µ sb − µ cb ) e − ( γ b /ν f )Ω √ GJ I which is typically positive due to dominance of the rock-bit over the mud-bit interaction.In contrast, the curvature parameter p = ψ (cid:48)(cid:48) (0) = − LI φ (cid:48)(cid:48) (Ω) = − W ob R b ( µ sb − µ cb ) ( γ b /ν f ) e − ( γ b /ν f )Ω LI is typically negative. 3. Analysis of the linear system G In this section we determine the number of unstable poles of the linearization G of G nl ,as this will be needed later to assure stability of the closed loop. This discussion is ofindependent interest, as in a different context the specific form of the non-linearity ψ ( x t ) may be unknown, in which case a linear parametric robust synthesis in q may be required.We recall from [18] that the open loop G nl is stable for φ (cid:48) > , and the same is truefor its linearization G . This means that we may concentrate on the potentially instablecase φ (cid:48) ≤ , which means q ≥ . Our goal is to classify the open loop properties of G asa function of the three parameters ( q, α, λ ) ∈ R .Laplace transformation of (1) leads to a family of one-dimensional boundary valueproblems parametrized by s ∈ C : G : x ξξ ( ξ, s ) = ( s + 2 λs ) x ( ξ, s ) x ξ (1 , s ) = − sx (1 , s ) + u ( s ) x ξ (0 , s ) = ( αs − qs ) x (0 , s ) (9) PIERRE APKARIAN AND DOMINIKUS NOLL Five scenarios of physical parametersgray blue magenta red green G N · m − J m I kg · m I B
89 35.6 35.6 35.6 89 kg · m L mΩ
10 10 10 10 10 rad · s − c a kg · m · s − β N · s W ob N R b m µ sb µ cb γ b ν f rad · s − c b Table 1.
Derived parameters of five scenariosgray blue magenta red green Ω rad · s − − Ω L √ I √ GJ -3.7186 -6.5044 -3.7186 -6.5044 -3.7186 rad λ α q p -0.0048 -0.1506 -0.2927 -0.1931 -0.2927 – τ /t s − n p n z Table 2. which we solve explicitly. With the outputs y ( s ) = sx (0 , s ) , y ( s ) = sx (1 , s ) from (2) weobtain G ( s ) = y ( s ) u ( s ) y ( s ) u ( s ) = e σ − e − σ σ (cid:2) σ s + αs − qs (cid:3) + e σ + e − σ [ αs − q + 1] e σ + e − σ + ( αs − qs ) e σ − e − σ σe σ − e − σ σ (cid:2) σ s + αs − qs (cid:3) + e σ + e − σ [ αs − q + 1] , σ ( s ) := √ s + 2 λs. (10)We now have to determine the number of unstable poles of (10) as a function of ( q, α, λ ) ∈ R . Note that G ( s ) = [1 /d ( s ); n ( s ) /d ( s )] is a meromorphic function, with n ( s ) , d ( s ) in(10) holomorphic, but its analysis is more complicated than that of a pure delay systemdue to the damping coefficient λ and the consequent appearance of the term σ ( s ) . ONTROL OF AN ANTI-STABLE WAVE 7
Annihilating d ( s ) = ( s + 2 λ + αs − qs ) e σ − e − σ σ + ( αs − q + 1) e σ + e − σ = 0 leads to thecomplex equation(11) q − αs − λs + ( e σ + e − σ ) / e σ − e − σ ) / σ =: Φ( λ, s ) , which relates unstable pole s ∈ C + of G and damping coefficient λ > to the pair ( q, α ) through the operator Φ . Since this is a complex equation and q, α are real, we deduce(12) α = − Im Φ( λ, s )Im( s ) , q − λ, s )Im( s ) − Im Φ( λ, s )Re( s )Im( s ) . We have proved the following
Lemma 1.
Let λ > and s ∈ C + . Suppose ( q, α ) given by (12) is in R . Then s is anunstable pole of G for the parameters ( q, α, λ ) ∈ R . (cid:3) Let us look at poles on the imaginary axis j R , referred to as zero-crossings . Going backto (11) with s = jω gives Lemma 2.
Let λ > and ω ∈ R . Suppose the pair q − λ, jω ) , α = − Im Φ( λ, jω ) ω satisfies ( q, α ) ∈ R . Then jω is a zero crossing (unstable pole on j R ) of G for theparameter ( q, α, λ ) ∈ R . (cid:3) Let us look more specifically at zero-crossings through the origin. Substituting s = 0 in the denominator d ( s ) in (10) and equating d (0) = 0 gives the relation q − λ crit which says that for q > a real pole of G crosses the imaginary axis through the originat the critical value λ = λ crit . Here we use the fact that Φ( λ,
0) = 2 λ , explained by therelations(13) e σ − e − σ σ = 1 + σ + σ + . . . , e σ + e − σ = 1 + σ + σ + . . . , σ s = s + 2 λ. Theorem 1.
For fixed ( q, α, λ ) ∈ R there exists R > such that G has no poles and notransmission zeros on { s ∈ C + : | s | ≥ R } . Proof:
1) For unstable poles we have to show that equation (11), respectively, (12) hasno solutions when s ∈ C + and | s | (cid:29) sufficiently large. Let s = µ + jω , σ = a + jb , thenby the definition of σ :(14) a − b = µ − ω + 2 λµ, ab = ω ( µ + λ ) . It follows that for fixed a > the set { s ∈ C + : Re( σ ) ≤ a } is bounded. Choose R > such that { s ∈ C + : Re( σ ) ≤ a } ⊂ { s ∈ C + : | s | ≤ R } . It remains to discuss candidatepoles s ∈ C + with Re ( σ ) ≥ a for some fixed a > .2) Consider s ∈ C + with Re ( σ ) = a ≥ a and define θ := e σ + e − σ e σ − e − σ = 1 + e − a e − j b − e − a e − j b . Then ρ − ρ ≥ | θ | ≥ − ρ ρ , where ρ = e − a . Moreover, we have | θ + 1 | ≥ e − a =: θ > .Now choose (cid:15) > such that ρ − ρ (cid:15) < θ / . Since λ is fixed we have σ/s → as s → ∞ on C + , hence there exists M = M ( λ ) > such that | s − σ | < (cid:15) | s | for all | s | ≥ M . Then PIERRE APKARIAN AND DOMINIKUS NOLL | s + θσ | = | θ ( σ − s ) + ( θ + 1) s | ≥ | θ + 1 || s | − | θ || s − σ | ≥ θ | s | − | θ | (cid:15) | s | ≥ θ / | s | for | s | ≥ M . Writing (11) as(15) q − − αs = Φ( λ, s ) = 2 λs + θσ , and taking into account that on the right hand side we now have | Φ( λ, s ) | = 2 λ | s + θσ | ≤ λθ | s | we see that (15) can have no solution for | s | ≥ max { λαθ , qα , M, R } =: R . Thatsettles the case α > .3) For α = 0 and q (cid:54) = 1 there are no poles in | s | > λ | q − | θ , and for q = 1 , α = 0 clearly(15) has no solutions.4) Let us next discuss unstable zeros. Clearly those can only occur in the secondcomponent y /u in (10). Here the equation is Φ( λ, s ) − = ( q +1) s − αs λ . From 1) above weknow that we may concentrate on Re ( σ ) ≥ a , and from 2) we have | θ | ≤ ρ − ρ , while for | s | > λ we get | σ | ≤ √ | s | , hence for | s | > max { R , λ } : | Φ( λ, s ) | = (cid:12)(cid:12)(cid:12)(cid:12) λs + θσ (cid:12)(cid:12)(cid:12)(cid:12) ≥ λ | s | + | θ || σ | ≥ λ | s | (1 + ρ − ρ √ . This leads to | ( q + 1) s − αs | λ = | Φ( λ, s ) | − ≤ | s | (1 + ρ − ρ √ λ hence α | s | − ( q + 1) ≤ | αs − ( q + 1) | ≤ ρ − ρ √ . For α > this cannot be satisfied for large | s | . In fact, there are no unstable zeros on | s | > R := max { R , λ, (2 + q + ρ − ρ √ /α } . For α = 0 the equation for unstable zerosis θσ = qs , and since σ/s = (cid:112) λ/s → for s → ∞ , we get θ → q . On choosing a sufficiently large, we get θ ≈ , which leads to a contradiction for q (cid:54) = 1 . Finally, for q = 1 , α = 0 we obtain the transfer function y /u = s [( σ − s ) e σ +( σ + s ) e − σ ]( σ − s )( σ + s ) e σ − ( σ + s ) e − σ , so unstablezeros (cid:54) = 0 satisfy s − σs + σ = e − σ . That gives − λs + σ + λ = e − σ , which cannot be satisfied forlarge | s | . (cid:3) Since the transfer function G is of size × , the number of unstable poles is themaximum of the number of unstable poles of G ( s ) = 1 /d ( s ) and G ( s ) = n ( s ) /d ( s ) ,hence the number of unstable zeros of d ( s ) . The latter can be determined by the argumentprinciple. For the following we denote the half circle used for the standard Nyquist contourby D R . Proposition 1.
Suppose ( q, α, λ ) ∈ R does not give rise to zero crossings. Then thenumber n p of unstable poles of G ( s ) equals the winding number of d ( D R ) around , wherethe radius R > is as in Theorem . The radius R in Theorem 1 may be quantified, and the winding number can be com-puted exactly using the method in [5]. If ( q, α, λ ) creates a zero-crossing, the contour D R has to be modified, either by making small indentations into the right half plane, orpreferably by removing poles on j R with the method of [25], as explained in [5]. At thisstage we have completed step 1 of our general algorithm 1.We conclude this section with the following important consequence of Theorem 1. ONTROL OF AN ANTI-STABLE WAVE 9
Corollary 1.
The input-output map and the input-to-state map of the boundary controlproblem (1) are bounded.
Proof:
As a consequence of [20, Thm. 2.3] for input-output boundedness it sufficesto show that sup
Re( s ) >a | G ( s ) | < ∞ for some a ∈ R . We choose a as in the proof ofTheorem 1, which allows to bring θ as close to 1 as we wish. Now with the notation ofthe theorem G ( s ) = σθ + αs − qsσ /s + αs − qs + 2 σθ ( αs − q + 1) . For α > we divide numerator and denominator by the leading term αs , which gives G ( s ) = 1 + σθ/αs − q/αs /αs + 2 λ/αs − q/αs + 2 θσ/s + 2(1 − q ) σθ/αs ∼
11 + 2 θσ/s .
But σ /s = 1 + 2 λ/s ∼ , whence G ( s ) ∼ / , showing that G is bounded on somehalf plane Re ( s ) > a . Since G = n/d and G = 1 /d , this is also true for G . In the case α = 0 simplification by s leads to a similar estimate. (cid:3) Pattern of unstable poles
As a consequence of the previous section we can determine the number n p of unstablepoles of G for every scenario ( q, α, λ ) ∈ R using the argument principle. However,we would like to learn a little more about n p ( q, α, λ ) , and in this section we shall seethat n p ∈ { , , } , where the corresponding regions can be determined with arbitrarynumerical precision.To begin with, observe that for λ = 0 the transfer function (due to σ = s ) simplifies toa pure delay system G λ =0 ( s ) = e − s αs − q (1 + αs − q ) + (1 − αs + q ) e − s αs − q ) = α e − s s − q − α
12 + (cid:0) qα − s (cid:1) e − s s − q − α , where we immediately see that G λ =0 has one unstable real pole if q ≥ , while it is stablefor q < .This suggests now the following procedure. Fix ( q, α ) ∈ R , and then follow theevolution of the number of unstable poles n p ( λ ) := n p ( q, α, λ ) of G as λ increases from λ = 0 to λ → + ∞ . We know the number of poles at λ = 0 , and we expect that forvery large λ (cid:29) the damping effect in the wave equation should lead back to stability, n p ( λ (cid:29)
0) = 0 .Let us look again at zero crossings at the origin. We know that for q > the originis crossed when λ ∈ [0 , ∞ ) reaches the critical value λ crit = ( q − / > . We have todecide whether this real pole when crossing s = 0 migrates from left to right or in theopposite direction. Let s ( λ ) be the position of the potentially unstable pole on the realaxis, that is d ( s ( λ ) , λ ) = 0 , where s ( λ crit ) = 0 . Differentiation with respect to λ gives s (cid:48) ( λ ) = − d λ ( s ( λ ) , λ ) d s ( s ( λ ) , λ ) , where ∂d∂λ = 2 (cid:18) σ
3! + . . . (cid:19) +( s +2 λ + αs − qs ) sσ (cid:18) σ
3! + 4 σ
5! + . . . (cid:19) +( αs − q +1) sσ e σ + e − σ AND DOMINIKUS NOLL and ∂d∂s = (1 + 2 αs − q ) (cid:18) σ
3! + . . . (cid:19) + ( s + 2 λ + αs − qs ) s + λσ (cid:18) σ
3! + 4 σ
5! + . . . (cid:19) + α e σ + e − σ αs − q + 1) e σ − e − σ σ ( s + λ ) . Substituting λ = λ crit = ( q − / and s = s ( λ crit ) = 0 gives s (cid:48) ( λ crit ) = 2 ( q − + ( q − − α . Hence s (cid:48) ( λ crit ) (cid:26) > for α < ( q − + ( q − < for α > ( q − + ( q − This leads to the following
Lemma 3.
Let ( q, α ) ∈ R . If α < ( q − + ( q − , then a single real pole of G crossesthe imaginary axis through the origin at λ = λ crit = ( q − / from left to right, goingfrom stable at λ crit − to unstable at λ crit + 0 . If α > ( q − + ( q − a single realpole crosses the imaginary axis through the origin from right to left, going from unstableat λ = λ crit − to stable at λ = λ crit + 0 . (cid:3) This can also be corroborated by investigating the value G (0) in (10). We have G (0) = 12 λ − q + 1 , so G has no unstable pole at the origin, except for the critical λ value λ crit = ( q − / when q > . On the exceptional manifold M = { ( q, α, λ ) ∈ R : 2 λ = q − } , we have lim s → sG ( s ) = 1 α − ( q − − ( q − , which means a pole of order one at the origin, except when ( q, α ) lies on the parabola α = ( q −
1) + ( q − . On the exceptional set O = { ( q, α, λ ) ∈ R : 2 λ = q − , α =( q −
1) + ( q − } we find that lim s → s G ( s ) = 6 q − , which means G has a double pole at the origin, except when q = 1 . The case q = 1 now leaves only the parameter choice ( q, α, λ ) = (1 , , , an exceptional point where thesystem is not well-posed.Using the mapping Φ , one can see that the positive quadrant ( q, α ) ∈ R may bedivided into 5 different zones, shown in Fig. 3, in which the number of unstable poles of G evolves differently. Each zone has its specific pattern.The red zone is Red = { ( q, α ) : q ≥ , α ≥ , α ≤ ( q − + ( q − } is below aparabola. Setting m ( q ) = sup { α > q − λ, jω ) , α = − ω − Im Φ( λ, jω ) for certain ω > , λ > } , the magenta zone is defined as Mag = { ( α, q ) : q ≥ , ( q − + ( q − ≤ α ≤ m ( q ) } delimited by the parabola and the analytic curve α = m ( q ) . The green zone is Green = { ( q, α ) : q ≥ , α ≥ m ( q ) } , ONTROL OF AN ANTI-STABLE WAVE 11 where the curve α = m ( q ) separates magenta and green. Finally, on setting b ( α ) = inf { q : q − λ, jω ) , α = − ω − Im Φ( λ, jω ) for certain ω > , λ > } , the blue zone is Blue = { ( q, α ) : α ≥ , b ( α ) ≤ q ≤ } , which is the only bounded one.The boundary of the blue zone described by the curve q = b ( α ) is just a different localparametrization of the same analytic curve α = m ( q ) separating magenta and green. Thiscurve disappears into α < at (1 , , where it is no longer of interest. The gray zone Gray is what is left over from the strip ≤ q ≤ , α ≥ when removing the blue zone. Ω angular velocity f r i c t i o n f o r ce φ Ω angular velocity f r i c t i o n f o r ce φ Figure 1.
The case (⌦) > (left) leads to a stable open loop [ ? ]. Thepotentially unstable case (right) is when increasing rotary speed reducesfriction. q α Figure 2.
Five regions and five scenarios. The pattern of the gray zone is0. Blue: 0-2-0. Red: 1-2-0. Magenta: 1-0-2-0. Green: 1-0. ψ ( · ) G y K u y W u z u W z Figure 3.
Five regions and five scenarios. The pattern of the gray zone is0. Blue: 0-2-0. Red: 1-2-0. Magenta: 1-0-2-0. Green: 1-0.Altogether, we have found the following classification or pattern. • For ( q, α ) ∈ Gray the system G is stable for all λ ≥ . The pattern is . • For ( q, α ) ∈ Blue there exist < λ ( q, α ) < λ ( q, α ) such that G is stable for all ≤ λ < λ ( q, α ) and λ > λ ( q, α ) , and has two unstable poles for λ ≤ λ ≤ λ .The pattern is . • For ( q, α ) ∈ Red the system has one unstable pole for ≤ λ ≤ ( q − / λ ( q ) ,and two unstable poles for ( q − / ≤ λ ≤ λ ( q, α ) , while it is again stable for λ > λ ( q, α ) . The pattern is . • For ( q, α ) ∈ Mag there exist λ ( q, α ) > λ ( q, α ) > ( q − / such that the systemhas one unstable pole for ≤ λ ≤ ( q − / , no unstable poles for ( q − / <λ < λ ( q, α ) , then two unstable poles for λ ( q, α ) ≤ λ ≤ λ ( q, α ) , and again nounstable poles for λ > λ ( q, α ) . The pattern is . • For ( q, α ) ∈ Green the system has one unstable pole for ≤ λ ≤ ( q − / , and isstable for λ > ( q − / . The pattern is .5. Stabilization
In this section we construct finite-dimensional output feedback controllers which stabi-lize the linearization G of system (1)–(2) exponentially. We start with the following Theorem 2.
Let K be a finite-dimensional output feedback controller for (1)-(2) whichstabilizes the closed loop in the H ∞ -sense. Then the closed loop is even exponentiallystable. AND DOMINIKUS NOLL Proof:
1) Suppose the boundary control problem (BCP) is written in the abstract form ˙ x = A x, P x = u, y = C x with suitable unbounded operators [20, 36, 37], and let the controller u = Ky stabilizeBCP in the H ∞ -sense. Writing K ( s ) = K ( s ) + K with K strictly proper, we see that (cid:101) u = K y stabilizes the modified BCP ˙ x = A x, ( P − K C ) x = (cid:101) u, y = C x in the H ∞ -sense, where (cid:101) u = u − K y . We will use this type of shift to arrange for a strictlyproper stabilizing controller.2) We start from (1) by performing the change of variables z ( ξ, t ) = x ξ ( ξ, t ) , v ( t ) = x t (0 , t ) , cf. [35], which leads to an equivalent representation of (1) as a PDE coupled withand ODE: G : z tt ( ξ, t ) = z ξξ ( ξ, t ) − λz t ( ξ, t ) z (1 , t ) = (cid:101) u ( t ) αz ξ (0 , t ) = z (0 , t ) + ( q + 2 αλ ) v ( t ) α ˙ v ( t ) = z (0 , t ) + qv ( t ) (16)where the new state is ( z, z t , v ) , the measured outputs are y = v, y ( t ) = (cid:90) z t ( ξ, t ) dξ + v ( t ) , and where a new control (cid:101) u ( t ) = u ( t ) − x t (1 , t ) = u ( t ) − y ( t ) is used. Since the controller u = Ky stabilizes (1) in the H ∞ sense by hypothesis, so does (cid:101) u = Ky − y for (16), andsince the state trajectories remain unaffected, we may from here on prove the statementfor controller (cid:101) u = Ky and system (16). It is also clear that we may replace the outputs y , y by equivalent outputs (cid:101) y = v , (cid:101) y = (cid:82) z t ( ξ, t ) dξ , because (cid:101) y = y , (cid:101) y = y − y .Then (cid:101) u = u − y = u + (cid:101) y − (cid:101) y , and the controller is (cid:101) u = (cid:101) K (cid:101) y . At this stage for the easeof presentation we drop the tilde notation and write the new control and measurementsagain as u and y .3) Let the controller K have the form u ( s ) = K ( s ) y ( s ) = K ( s ) y + K y with directtransmission K y = k y + k y and strictly proper part K ( s ) . We now apply the ideaof part 1) and shift its direct transmission into the plant. This leads to G (cid:48) : z tt − z ξξ + 2 λz t = 0 z (1 , t ) = u − k y − k y αz ξ (0 , t ) − z (0 , t ) = ( q + 2 αλ ) v ( t ) α ˙ v = qv + z (0 , t ) (17)with the outputs as before, now in feedback with u = K ( s ) y . Note that K still stabilizes(17) in the H ∞ -sense, and since the state trajectories remain the same we may proveexponential stability of the loop for this pair G (cid:48) , K . Since K is strictly proper, thecontroller ˙ u = K y , respectively, u = s K ( s ) y = K (cid:48) ( s ) y , is proper and may be representedin state space as K (cid:48) : ˙ x K = A K x K + B K y + B K y ˙ u = (cid:101) u (cid:101) u = C K x K + d K y + d K y . ONTROL OF AN ANTI-STABLE WAVE 13
If the original state-space realization is K = (cid:20) a bc d (cid:21) , then K = (cid:20) a bc (cid:21) , and K (cid:48) = (cid:20) a bca cb (cid:21) =: (cid:20) A K B K C K D K (cid:21) . Since H ∞ -stability of the loop is not altered by these trans-formations, we may prove the statement for the pair G (cid:48) , K (cid:48) .4) We now perform a less standard manipulation, which consists in transferring partsof the system dynamics (17) into a new augmented controller (cid:101) K . We introduce a newartificial output y = z (0 , t ) in (17), and consider the boundary wave equation (cid:101) G : z tt − z ξξ + 2 λz t = 0 z (1 , t ) = u ( t ) − k v − k y αz ξ (0 , t ) − z (0 , t ) = ( q + 2 αλ ) v ( t ) y ( t ) = (cid:90) z t ( ξ, t ) dξy ( t ) = z (0 , t ) (18)Here we have substituted v = y , created a new input into (cid:101) G , and have now an infinitedimensional system (cid:101) G in feedback with the extended controller (cid:101) K : α ˙ v = qv + y ˙ x K = A K x K + B K v + B K y ˙ u = C K x K + d K v + d K y ˙ v = qα v + 1 α y (19)The ODE α ˙ v = z (0 , t ) + qv = qv + y was shifted from G (cid:48) into the new (cid:101) K , leaving us witha simpler infinite-dimensional system (cid:101) G . The controller (cid:101) K is K (cid:48) augmented by this ODE,so is still finite dimensional, and moreover, is also an integral controller with regard toits new output v . The output y has disappeared from (18), because the correspondingdynamics are now integrated in (cid:101) K . The state of (18) is ( z, z t ) , while the state of (cid:101) K is ( v, x K ) , to which we have to add the integrator. (cid:101) K is an integral controller with regardto the new output ˙ v from (19).5) Our next step is to find a state-space representation of y = (cid:101) G [ u ; v ] T in (18), whichmeans representing it as a well-posed boundary control system in the sense of [36, 37],[44, Def. 5.2.1] or [20]. With zero boundary conditions equation (18) reads z tt − z ξξ + 2 λz t = 0 z (1 , t ) + k (cid:90) z t ( ξ, t ) dξ = 0 αz ξ (0 , t ) − z (0 , t ) = 0 z ( ξ,
0) = z ( ξ ) , z t ( ξ,
0) = z ( t ) . This has now a representation as a strongly continuous semi-group(20) ˙ z = (cid:20) I d dξ − λ (cid:21) z =: A z , z (0) = z , where z = ( z, z t ) , and where the generator A has D ( A ) = { ( z , z ) ∈ H (0 , × H (0 ,
1) : z (1) + k (cid:82) z ( ξ, t ) dξ = 0 , αz x (0) − z (0) = 0 } as domain in the Hilbert space H = H (0 , × L (0 , . Define A with domain D ( A ) = H × H by the same formula (20), and AND DOMINIKUS NOLL let the projector P with D ( P ) = D ( A ) be defined as P z = (cid:20) z (1) + k (cid:82) z t ( ξ ) dξαz ξ (0) − z (0) (cid:21) ∈ C .The boundary control of (cid:101) G has now the abstract form ˙ z = A z , P z = (cid:20) u − k vk v (cid:21) , y = C z , where k := q + 2 αλ , and where y = [ y , y ] , with C : H × L → C bounded. Finallywe re-arrange the boundary condition by defining P (cid:48) with D ( P (cid:48) ) = D ( P ) as P (cid:48) z = (cid:34) z (1) + k (cid:82) z t ( ξ ) dξ + α k k z ξ (0) − k k z (0) k ( αz ξ (0) − z (0)) (cid:35) , P (cid:48) z = (cid:20) uv (cid:21) =: u . In order to make this well defined, we have according to [22, Sect. 3.3] to assure that D ( A ) ⊂ D ( P (cid:48) ) , D ( A ) = D ( A ) ∩ ker( P (cid:48) ) , A z = A z on D ( A ) , and that A generates a C -semi group. These are satisfied by construction. In addition, we require an operator B : C → H such that P (cid:48) ◦ B = I , im ( B ) ⊂ D ( A ) and A ◦ B bounded. This can bedefined by the ansatz B u = (cid:20) b ( ξ ) u + c ( ξ ) v (cid:21) , b ( ξ ) = ξ , c ( ξ ) = ( k − k ) ξ − k . Indeed, P (cid:48) B u = (cid:20) b (1) u + c (1) v + d K ( αb (cid:48) (0) u + αc (cid:48) (0) v − b (0) u − c (0) v ) αb (cid:48) (0) u + αc (cid:48) (0) v − b (0) u − c (0) v (cid:21) = (cid:20) uv (cid:21) = u . Asis well-known, (cf. [22, Sect. 3.3]), the boundary wave equation may now be representedby the state-space(21) ˙ x = A x − B ˙ u + A B u , x (0) = x , where solutions z of (18) and x of (21) are related by x ( t ) = z ( t ) − B u ( t ) . This can befurther streamlined as(22) ˙ x e = (cid:20) A B A (cid:21) x e + (cid:20) I − B (cid:21) (cid:101) u , where the extended state is x e = ( u , x ) = ( u, v, z, z t ) , and where (cid:101) u = ˙ u has become theinput. The output operator for (22) is now y = C e x e = C ◦ [ B I ] (cid:20) u z − Bu (cid:21) = C z .6) We next show that system (cid:101) G , and therefore also the state-space representation (22)with C -semi group, is exponentially stabilizable. This can for instance be obtained from[35, Theorem 4.2], where the authors construct a state feedback controller which stabilizes(16) exponentially in the Hilbert space H = H × H × L . The control law found in thatreference can be arranged as a state feedback law for (cid:101) G , and hence for (22), using thesame technique of shifting parts of the dynamics from plant to controller. Alternatively,we may even use the open loop characterization of stabilizability, called optimizability in[48], which is equivalent to stabilizability, while offering a more convenient way to checkit.7) We now show that the controller (cid:101) K is admissible for (cid:101) G and is as a system exponen-tially stabilizable. Due to shifting the direct transmission of K into the plant as outlinedin 1) and put to work in (18), the new controller (cid:101) K in (19) is written as an integralcontroller, that is, its output is (cid:101) u = ˙ u = [ ˙ u, ˙ v ] , which makes it admissible for (cid:101) G .Assuming that the original K = (cid:20) a bc d (cid:21) is stabilizable and detectable (e.g. mini-mal), the same is true for K (cid:48) obtained in 3), so (cid:20) A K B K C K D K (cid:21) is stabilizable. Now the ONTROL OF AN ANTI-STABLE WAVE 15 augmented controller is (cid:101) K = (cid:34) (cid:101) A K (cid:101) B K (cid:101) C K (cid:101) D K (cid:35) with (cid:101) A K = (cid:20) A K B K α − (cid:21) , (cid:101) B K = (cid:20) B K α − (cid:21) , (cid:101) C K = (cid:20) C K d K α − (cid:21) , (cid:101) D k = (cid:20) d K α − (cid:21) . Applying the Hautus test, for simplicity in the case λ (cid:54) = α − , let v be an eigenvector of A TK with unstable eigenvalue λ , then [ v ρ ] T is aneigenvector of (cid:101) A TK for λ if ρ = B TK v/ ( λ − α − ) . Now (cid:101) B TK [ v ρ ] T = [ v v ( α − / ( λ − α − ))] T ,and this vector cannot be = [0 0] T , because that would imply B TK v = 0 and B TK v = 0 ,hence B TK v = 0 , contradicting stabilizability of [ A K , B K , C K , D K ] . Now for the eigenvalue α − of (cid:101) A TK we take the eigenvector w = [0 1] T , then (cid:101) B TK w = [0 α − ] T (cid:54) = [0 0] T , whichproves stabilizability.With (cid:101) G and (cid:101) K exponentially stabilizable, the closed loop ( (cid:101) G, (cid:101) K ) is also exponentiallystabilizable (see [44, Prop. 8.2.10(ii)(c)]) in the sense of the induced state-space realization[44, Chap. 7]. The infinitesimal generator of the closed loop will be denoted as A cl .8) Next we argue that (cid:101) G is exponentially detectable. Since (cid:101) G is exponentially stabiliz-able, its semi-group satisfies the spectrum decomposition assumption, see [22, Theorem5.2.6]. Since from the discussion of section 3 we know that there are only finitely manyright hand poles, all with finite multiplicity, a necessary and sufficient condition for ex-ponential detectability is that ker ( sI − (cid:101) A ) ∩ ker( (cid:101) C ) = { } for every s ∈ C + ; see [22,Theorem 5.2.11], where ( (cid:101) A, (cid:101) B, (cid:101) C ) refers to the state-space realization of (cid:101) G derived in 2)above. But that may now be checked in the frequency domain. It means that for every s ∈ C + the only solution of the Laplace transformed system (9) with u = 0 satisfying y ( s ) = sx (0 , s ) = 0 , y ( s ) = sx (1 , s ) = 0 is x ≡ . Now for s (cid:54) = 0 these boundaryconditions give x (0 , s ) = 0 , x (1 , s ) = 0 , and therefore from the boundary conditions in(9) x ξ (0 , s ) = 0 , x ξ (1 , s ) = 0 . The general solution of the dynamic equation in (9) be-ing x ( ξ, s ) = k e σξ + k e − σξ , with constants depending on s , we get the four conditions k + k = 0 , σ ( k − k ) = 0 , k e σ + k e − σ = 0 , σ ( k e σ − k e − σ ) = 0 , which can only besatisfied if k = k = 0 .With (cid:101) G exponentially detectable, and (cid:101) K exponentially detectable with an argumentsimilar to 7) above, the closed loop is exponentially detectable, again by [44, Prop.8.2.10(ii)(c)].9) According to [26, Theorem 5.2], a well-posed system which is exponentially stabiliz-able, exponentially detectable, and at the same time H ∞ -stable, is already exponentiallystable in the state-space sense, i.e., the generator of its semi-group is exponentially stable.We apply this to the closed loop system with generator A cl . For this result see also [48,Thm. 1.1], and [29, 8.35] for a classical antecedent. (cid:3) Corollary 2.
Let K be a finite-dimensional controller for (1) - (2) , and suppose the closedloop with G has no unstable poles. Then K stabilizes G exponentially and G nl locallyexponentially. Proof:
The result follows from Theorem 2 above once we show that (cid:101) K in (19) stabi-lizes (cid:101) G in (16) in the H ∞ -sense. Since the transfer functions are not concerned by thetransformations in the proof of theorem 2, it suffices to show that K stabilizes G in the H ∞ -sense. For that we have to show that the closed loop transfer operator T ( s ) = (cid:20) I G ( s ) − K ( s ) I (cid:21) − = (cid:20) ( I + KG ) − − K ( I + GK ) − ( I + GK ) − G ( I + GK ) − (cid:21) ∈ H ∞ belongs to the Hardy space H ∞ . Since we know by hypothesis that T ( s ) has no polesin C + , this follows as soon as T is bounded on j R . Since K is proper, this hinges AND DOMINIKUS NOLL on the behavior of G on j R . As is easy to see, the denominator d ( s ) of (10) satisfies lim ω →∞ | d ( jω ) | = ∞ , so y ( s ) /u ( s ) is proper, and it remains to show that y ( s ) /u ( s ) in (10) is bounded on j R . Dividing numerator and denominator of n /d in (9) by ( e σ − e − σ ) / σ , and observing that the term ( e σ + e − σ ) / e σ − e − σ ) / σ is bounded on j R , we see that y ( s ) /u ( s ) is bounded, because the leading terms in both numerator and denominator are now αs = − αω . That proves H ∞ -stability of the closed loop hence exponential stability of the linearclosed loop.It remains to show that K stabilizes G nl locally exponentially. Due to the specific formof the non-linearity, this may be obtained with [51]. (cid:3) Remark 1.
Semi-groups for hyperbolic equations with boundary dynamics have beeninvestigated, e.g. in [31], but as this requires additional conditions, we believe that ourmethod of simplifying the infinite-dimensional part by augmenting the controller offersadditional flexibility. After all the goal is to show that the closed-loop has an exponentiallystabilizable and detectable semi-group e A cl t , not necessarily the individual parts.We have achieved that for any finite-dimensional controller K for (1)-(2), exponentialstability of the infinite dimensional closed loop with G can be verified by the Nyquist test.What remains to do is actually find such a stabilizing controller. A straightforward ideais to use a discretization of (1), the most obvious being finite differences x i ( t ) = x ( ξ i , t ) , ξ i = ih, i = 0 , . . . , N, N h = 1 x ξ ( ξ i , t ) ≈ x i +1 ( t ) − x i − ( t )2 h , x ξξ ( ξ i , t ) ≈ x i +1 ( t ) + x i − ( t ) − x i ( t ) h . With the boundary condition at ξ = 0 αx (cid:48)(cid:48) ( t ) = x ( t ) − x − ( t )2 h + qx (cid:48) ( t ) we can eliminate x − , and with the boundary condition at ξ = 1 x N +1 ( t ) − x N − ( t )2 h = − x (cid:48) N ( t ) + u ( t ) we eliminate x N +1 . Putting ˜ x i = x (cid:48) i , i = 0 , . . . , N , we get a dynamical system of order N + 2 (23) (cid:20) x (cid:48) ˜ x (cid:48) (cid:21) = (cid:20) IT Λ (cid:21) (cid:20) x ˜ x (cid:21) + (cid:20) b (cid:21) u, y ( t ) = ˜ x ( t ) , y ( t ) = ˜ x N ( t ) , with typical A -matrix featuring a tridiagonal T and a diagonal Λ . It comes as a mildsurprise that (23) is not stabilizable, the reason being a pole/zero cancellation at theorigin.Kalman reduction using the function minreal from [52] removes one state of (23) andfurnishes a stabilizable system, which we use for synthesis, and where the reduced system A -matrix is now no longer sparse. In our experiment we chose N = 50 and synthesizedcontrollers of various simple structures like a sum of PIDs u = PID y + PID y , a 5th-order state-space controllers, or on ignoring one of the outputs, standard PID controllers u = PID y , respectively, u = PID y . These controllers, once they stabilize the reducedfinite-dimensional system, are then tested against the infinite dimensional system usingthe Nyquist test of [5], which gives an exact answer. For instance, stabilizing the grayand blue scenarios with the 5th-order controller given in (30) leads to the Nyquist plotsin Fig. 4 for the gray and blue scenarios, and certifies infinite-dimensional stability.As can be seen, in the blue case (right) the Nyquist curve winds twice around theorigin. Since K blue is stable and the open loop G blue has n p = 2 unstable poles, this proves ONTROL OF AN ANTI-STABLE WAVE 17
Gray scenario: Nyquist locus -25 -20 -15 -10 -5 0 5-15-10-5051015
Blue scenario: Nyquist locus
Figure 4.
Nyquist curve K blue G blue (right) winds twice around origin.Since K blue is stable and n p = 2 , closed loop is certified exponentially stable.Gray case (left) has n p = 0 and winding number around origin (criticalpoint). Since K gray is stable, loop is certified exponentially stable.exponential stability of the closed-loop ( G blue , K blue ) . At this point we have completedstep 2 of the general synthesis algorithm 1. Remark 2.
The fact that finite-difference and finite-element discretizations of stabilizable(or detectable) hyperbolic equations may turn out not stabilizable (detectable) cannotbe overcome by increasing N . This has been the cause of a large body of controllabilityliterature, which fortunately has little relevance for control. Namely, once we have decidedthat the true model for the drilling process is the infinite-dimensional (1)-(2), we littlecare whether K , synthesized for G and G nl , also stabilizes discretizations of G or G nl .6. H ∞ -synthesis The final step in algorithm 1 is H ∞ -synthesis. While we have already shown that thenon-linear system can be locally exponentially stabilized by a finite-dimensional controller,we now strive to prove global exponential stability of the closed-loop system ( G nl , K ) . Inthe following, it is helpful to represent the non-linear system G nl in Lur’e form, i.e., asthe closed loop interconnection of its linearization with a static non-linearity. K u G G nl ψ ( · ) y y w - - K u G y y w W u W y z y z u - - Figure 5.
Non-linear system (left) in feedback form. The synthesis inter-connection (right) interprets non-linearity as an exogenous disturbance w . AND DOMINIKUS NOLL Mixed sensitivity.
The non-linear system can be consider as a feedback loop be-tween the linearized plant P : x tt ( ξ, t ) = x ξξ ( ξ, t ) − λx t ( ξ, t ) x ξ (1 , t ) = − x t (1 , t ) + u ( t ) αx tt (0 , t ) = x ξ (0 , t ) + qx t (0 , t ) + w ( t ) y ( t ) = x t (0 , t ) , y ( t ) = x t (1 , t ) , z = ( y, u ) , (24)connected with the controller u = Ky and the non-linearity ψ ( · ) as in Fig. 5 left. Wenow have several choices. The most straightforward one is to grossly interpret the non-linearity ψ ( x t (0 , t )) as a disturbance w , forgetting its specific form. In (24) we thenintroduce typical outputs like z y = W y y , z u = W u u , where the channel w → z y rejectsthe effect of the non-linearity on the low-frequency part of the measured output, while w → W u u = z u accounts for high frequency components of the control signal, so thatminimizing the H ∞ -norm of T wz ( K ) limits the degrading effects of the non-linearity whilemaintaining reasonable control authority. Here and for the following T ab ( K ) denotes aclosed-loop channel b → a in plant P . The closed loop of (24) with K from w to z isobtained as T zw ( K ) = diag( W u , W y ) T ( u,y ) ,w ( K ) as shown in Fig. 5 right.6.2. Sector non-linearity.
A more sophisticated approach uses the fact that the non-linearity ψ in (7) induced by φ = φ mud + φ rock as in (8) is sectorial. That is to say, thereexist q l ≤ q u such that q l ω ≤ ψ ( ω ) ≤ q u ω for all ω , i.e., (1) is an infinite dimensional Lur’esystem. For the scenarios gray and blue these sectors are shown in Fig. 6. -8 -6 -4 -2 0 2 4 6-30-20-10010203040 Gray scenario: ( ) with sector bounds -10 -5 0 5 10-100-50050100150
Blue scenario: ( ) with sector bounds
Figure 6.
Sector non-linearity q l ω ≤ ψ ( ω ) ≤ q u ω for gray scenario (left)and blue scenario (right). Lemma 4.
For ω → ±∞ the non-linearity ψ ( ω ) behaves asymptotically like a line − q s ω + a ± , where q s = c b √ GJ I + q, a + = LW ob R b GJ ( µ sb − µ cb ) e − γbνf Ω , a − = 2 LW ob R b µ cb GJ + LGJ ( µ cb − µ sb ) e − γbνf Ω . Proof:
Note that since we have transferred the steady state to the origin, the kinkof the friction term ψ occurs at x t = − L Ω √ I √ GJ =: − x t . For x t (cid:29) − x t we have ψ ( x t ) = − ( c b √ GJI + q ) x t + LW ob R b GJ ( µ sb − µ cb ) e − γbνf Ω (1 − e − γbνf L √ GJI x t ) ∼ − ( c b √ GJI + q ) x t + LW ob R b GJ ( µ sb − µ cb ) e − γbνf Ω = − q s x t + a + , and for x t < x t we get ψ ( x t ) = − ( c b √ GJI + q ) x t + 2 LW ob R b µ cb GJ + LGJ ( µ cb − µ sb ) e − γbνf Ω + LGJ (cid:16) µ cb − µ sb ) e − L √ GJI | x t | e γbνf Ω (cid:17) ∼ − ( c b √ GJI + q ) x t + 2 LW ob R b µ cb GJ + LGJ ( µ cb − µ sb ) e − γbνf Ω = − q s + a − . (cid:3) ONTROL OF AN ANTI-STABLE WAVE 19
Since both branches behave asymptotically like a line with slope(25) − q s := − c b √ GJ I − q = − W ob R b ( γ b /ν f )( µ sb − µ cb ) e − γbνf Ω √ GJ I , it is not hard to find slopes q l , q u with q l ω ≤ ψ ( ω ) ≤ q u ω . Those can be seen in Fig. 6 forthe gray and blue cases. We use the standard notation ψ ∈ sect ( q l , q u ) .In order to achieve stability of the non-linear closed loop, we now apply the techniqueof Zames [49], which requires that the linear system T y w ( K ) in feedback with the non-linearity ψ ( · ) as in Fig. 5 right satisfy the complementary sector constraint. To put this towork, we let c = ( q l + q u ) / and r = ( q u − q l ) / , and introduce the centered non-linearity χ ( w ) = ψ ( w ) − cw , which satisfies χ ∈ sect ( − r, r ) .The centered non-linearity χ ( w ) = ψ ( w ) − cw is now in feedback with the followingshifted plant: (cid:101) P : x tt ( ξ, t ) = x ξξ ( ξ, t ) − λx t ( ξ, t ) x ξ (1 , t ) = − x t (1 , t ) + u ( t ) αx tt (0 , t ) = x ξ (0 , t ) + ( q + c ) x t (0 , t ) + e ( t ) y ( t ) = x t (0 , t ) , y ( t ) = x t (1 , t ) , z ( t ) = x t (0 , t ) , (26)connected with(27) u = Ky, z χ = χ ( e χ ) , e = z χ + w, e χ = z + w χ . Closing the loop with regard to u = Ky leads to z = (cid:101) T ze ( K ) e , which is in loop with thenon-linearity z χ = χ ( e χ ) as in Fig. 7. Here and in the following channels derived formplant (cid:101) P will be denoted (cid:101) T wz ( K ) etc. Note that the sole difference between P and (cid:101) P isthat the parameter q is replaced by (cid:101) q = q + c . In particular, stabilization of (cid:101) P is obtainedas studied in section 5. Ultimately this means that K will have to stabilize the linearwave equation for two different values q, (cid:101) q , while α, λ remain fixed. ! T ze ( K ) χ ( · ) z χ wze χ ew χ Figure 7.
Closing the loop with u = Ky in (26) leaves an exponentiallystable linear system (cid:101) T ze ( K ) in feedback with the shifted static non-linearity z χ = χ ( e χ ) . Lemma 5.
Let ψ ∈ sect ( q l , q u ) and put c = ( q u + q l ) / , r = ( q u − q l ) / . Sup-pose the controller K has been tuned such that the closed loop ( (cid:101) P , K ) is H ∞ -stable with (cid:107) (cid:101) T ze ( K ) (cid:107) ∞ < r − . Then the non-linear closed-loop (1) with u = Ky is finite gaininput-output stable, i.e., there exists a constant M > such that in (26)-(27) we have (cid:107) x t (0 , · ) (cid:107) + (cid:107) ψ ( x t (0 , · )) (cid:107) ≤ M ( (cid:107) w χ (cid:107) + (cid:107) w (cid:107) ) for all inputs w, w χ ∈ L [0 , ∞ ) . AND DOMINIKUS NOLL Proof:
This follows from [49, Thm. 1]. If ψ ∈ sect ( q l , q u ) , then the centered non-linearity χ := ψ − cI satisfies χ ∈ sect ( − r, r ) , hence (cid:107) χ ( z χ ) (cid:107) ≤ r (cid:107) z χ (cid:107) has L -gain r in the sense of [49, Def. (3)]. This non-linearity is now in feedback with (cid:101) T ze ( K ) .Since by assumption K has been tuned such that (cid:107) (cid:101) T ze ( K ) (cid:107) ∞ < r − , this LTI-systemhas L -gain < r − , and the small gain theorem implies boundedness of the loop (26)-(27), i.e., there exists a constant M > such that (cid:107) e χ (cid:107) ≤ M ( (cid:107) w (cid:107) + (cid:107) w χ (cid:107) ) and (cid:107) e (cid:107) ≤ M ( (cid:107) w (cid:107) + (cid:107) w χ (cid:107) ) in Fig. 7. Since in closed loop the input e to (cid:101) T ze ( K ) rep-resents the non-linear term χ ( x t ) , we derive (cid:107) χ ( x t (0 , · )) (cid:107) ≤ M ( (cid:107) w χ (cid:107) + (cid:107) w (cid:107) ) . Stillfrom the small gain theorem we get (cid:107) e χ (cid:107) ≤ M ( (cid:107) w χ (cid:107) + (cid:107) w (cid:107) ) , and since in closedloop e χ represents the output x t (0 , · ) , we have (cid:107) x t (0 , · ) (cid:107) ≤ M ( (cid:107) w χ (cid:107) + (cid:107) w (cid:107) ) in closedloop. Finally, for ψ = χ + cI we get a similar estimate by combining the previous two: (cid:107) ψ ( x t (0 , · )) (cid:107) ≤ (cid:107) χ ( x t (0 , · )) (cid:107) + c (cid:107) x t (0 , · ) (cid:107) ≤ M (1 + c ) ( (cid:107) w χ (cid:107) + (cid:107) w (cid:107) ) . (cid:3) This has now the following consequence:
Proposition 2.
Let ψ ∈ sect ( q l , q u ) with c, r as above, and suppose the controller K hasbeen tuned such that the closed loop ( (cid:101) P , K ) is H ∞ -stable with (cid:107) (cid:101) T ze ( K ) (cid:107) ∞ < r − . Thenthe non-linear closed loop between (1) and u = Ky is input-to-state stable in the followingsense: If the input signal w ∈ L [0 , ∞ ) , then the state ( x, x t ) of the non-linear closed loopwith initial condition x cl (0) = x is in L ([0 , ∞ ) , H ) . Proof:
Write the non-linear closed loop in the abstract state-space H in theorem 2 as ˙ x cl = A cl x cl +Ψ( x cl )+ w , x cl (0) = x , where A cl is exponentially stable, Ψ( x cl ) = ψ ( x t (0 , · )) for the closed-loop x t (0 , t ) , and where w ( t ) is an input to the equation αx tt (0 , t ) = x ξ (0 , t ) + qx t (0 , t ) + ψ ( x t (0 , t )) . The linear feedback system ( (cid:101) P , K ) , respectively itschannel (cid:101) T ze ( K ) , is now ˙ x cl = A cl x cl + e , z = C cl x cl , in loop with the centered non-linearity χ ( · ) , and w is the lower right input in Fig. 7. To account for a non-zeroinitial condition x cl (0) = x we choose the top left input in Fig. 7 as w χ = C cl e A cl t x ,where C cl is the output operator of closed loop system ( (cid:101) P , K ) . Then e χ is the solu-tion of the Cauchy problem ˙ x cl = A cl x cl + e, x cl (0) = x . From the lemma we get (cid:107) Ψ( x cl ) (cid:107) = (cid:107) ψ ( x t ) (cid:107) ≤ M ( (cid:107) w χ (cid:107) + (cid:107) w (cid:107) ) ≤ M (cid:48) ( (cid:107) x cl (0) e − ω t (cid:107) + (cid:107) w (cid:107) ) , where − ω < is the growth rate of the exponentially stable generator A cl . In particular, if we put v ( t ) = Ψ( x cl ( t )) + w ( t ) , then (cid:107) v (cid:107) ≤ ( M (cid:48) + 1) ( | x cl (0) | + (cid:107) w (cid:107) ) , hence we may consider v ( t ) as a right hand side in L to the non-homogeneous Cauchy problem ˙ x cl = A cl x cl + v , x cl (0) = x . Since A cl is exponentially stable, the closed loop state is then also in L ; [29,Ch. VI,7.1a]. (cid:3) One wonders whether the state x cl ( t ) decays exponentially to 0 when this is the casefor the input w ( t ) . Suppose w ∈ L decays exponentially in the sense that w = e − at (cid:101) w for some a > and (cid:101) w ∈ L . In this case it seems plausible to work with the weighted L -norm (cid:57) w (cid:57) = (cid:107) e at w ( · ) (cid:107) = (cid:107) (cid:101) w (cid:107) . Proposition 3.
Suppose ψ ∈ sect ( q l , q u ) with c, r as above, and suppose K has beentuned such that ( (cid:101) P , K ) is H ∞ -stable, with (cid:107) (cid:101) T ze ( K ) (cid:107) ∞ < r − . There exists a > suchthat whenever the input w decays exponentially with rate at least as fast as a , i.e., w ( t ) = e − at (cid:101) w ( t ) for some (cid:101) w ∈ L [0 , ∞ ) , then the state x cl ( t ) of the non-linear closed loop inresponse to the input w decays exponentially with rate at least a . Proof:
1) Since the closed loop ( (cid:101) P , K ) is exponentially stable with − ω := ω ( A cl ) < and (cid:107) (cid:101) T ze ( K ) (cid:107) ∞ < r − , we may choose a small enough shift < a < ω such that ONTROL OF AN ANTI-STABLE WAVE 21 ( (cid:101) P ( · − a ) , K ( · − a )) is still exponentially stable and (cid:107) (cid:101) T ze ( K )( · − a ) (cid:107) ∞ < r − . Let (cid:57) · (cid:57) be the corresponding weighted L -norm as above.2) Let us observe that for the centered non-linearity χ ∈ sect ( − r, r ) implies (cid:57) χ ( w ) (cid:57) ≤ r (cid:57) w (cid:57) for all w = e at (cid:101) w . Namely, (cid:57) χ ( w ) (cid:57) = (cid:82) t e aτ | χ ( w ( τ )) | dτ ≤ (cid:82) t e aτ r | w ( τ ) | dτ = r (cid:57) w (cid:57) .3) Now we establish the complementary estimate for the LTI feedback system ( (cid:101) P , K ) and its channel (cid:101) T ze ( K ) with regard to the norm (cid:57) · (cid:57) . We have (cid:57) (cid:101) T ze ( K ) ∗ w (cid:57) = (cid:90) ∞ e at (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) t (cid:101) T ze ( K )( t − τ ) w ( τ ) dτ (cid:12)(cid:12)(cid:12)(cid:12) dt = (cid:90) ∞ (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) t (cid:101) T ze ( K )( t − τ ) e a ( t − τ ) w ( τ ) e aτ dτ (cid:12)(cid:12)(cid:12)(cid:12) dt = (cid:90) ∞ (cid:90) t (cid:12)(cid:12)(cid:12)(cid:16) (cid:101) T ze ( K ) · e at (cid:17) ( t − τ ) (cid:0) w · e at (cid:1) ( τ ) dτ (cid:12)(cid:12)(cid:12) dt = (cid:107) ( (cid:101) T ze ( K ) · e at ) ∗ (cid:101) w (cid:107) = (cid:107) (cid:101) T ze ( K )( s − a ) · (cid:101) w ( s ) (cid:107) ≤ (cid:107) (cid:101) T ze ( K )( · − a ) (cid:107) ∞ (cid:107) (cid:101) w (cid:107) = (cid:107) (cid:101) T ze ( K )( · − a ) (cid:107) ∞ (cid:57) w (cid:57) < r − (cid:57) w (cid:57) . This means we may apply the small gain argument with the norm (cid:57) · (cid:57) . The result isas before that (cid:57) x t (0 , · ) (cid:57) + (cid:57) ψ ( x t (0 , · )) (cid:57) ≤ M ( (cid:57) w χ (cid:57) + (cid:57) w (cid:57) ) for some M > andall inputs w = e − at (cid:101) w , w χ = e − at (cid:101) w χ with (cid:101) w, (cid:101) w χ ∈ L [0 , ∞ ) . That means the non-linearityin closed loop in response to the signal w = e − at (cid:101) w also decays at least as fast as e − at , sothat the right hand side v ( t ) = Ψ( x cl ( t )) + w ( t ) already used in the previous propositionis of the form v ( t ) = e − at (cid:101) v ( t ) for some (cid:101) v ∈ L .We also have to argue that w χ = C cl e A cl t x decays with rate a , which holds since a < − ω ( A ) . But now all we have to observe is that due to exponential stabilityof A cl in the non-homogeneous Cauchy problem ˙ x cl = A cl x cl + v the state decays ex-ponentially as soon as v decays exponentially. The mild solution in the semi-groupsense [29, p. 436] satisfies x cl ( t ) = e A cl t x cl (0) + (cid:82) t e A cl ( t − τ ) v ( τ ) dτ , hence | x cl ( t ) | ≤ M (cid:16) e − ω t + (cid:107) (cid:101) v (cid:107) (cid:82) t e − ω ( t − τ ) e − aτ dτ (cid:17) ≤ M (1 + (cid:107) (cid:101) v (cid:107) / ( ω − a )) e − at . (cid:3) This brings us now to our first optimization program, where we combine a mixed H ∞ performance and robustness requirement (Fig. 5 right) for the nominal plant with a sectorconstraint assuring global exponential stability of the non-linear closed loop (Fig. 5 left)when satisfied: minimize r (cid:107) (cid:101) T ze ( K ) (cid:107) ∞ subject to (cid:107) W u T uw ( K ) (cid:107) ∞ ≤ K ∈ K (28)Here K refers to a class of structured controllers, and optimization over K ∈ K can bedispensed with as soon as the objective attains a value < . As our experiments show,the sectorial approach works successfully for the gray scenario. Note that it is implicit in(28) that K has to stabilize P and (cid:101) P , which means stabilizing the wave equation for thetwo different values q and (cid:101) q = q + c with the same α, λ .6.3. Large magnitude sector constraint.
The limitation of the sector approach isobviously that if the primal sector sect ( q l , q u ) is large, it is difficult to tune K suchthat the closed loop system ( P, K ) is in the complementary sector. In the transformed AND DOMINIKUS NOLL metric, if the primal sector is large, then r is large, so r − is small and the constraint (cid:107) (cid:101) T ze ( K ) (cid:107) ∞ < r − in (28) is difficult to achieve – if at all. This fails indeed for the bluescenario, and Zames-Falb multipliers [50] do not help for the specific non-linearity ψ .However, the particular structure of the non-linearity in Lemma 4 suggests the followingdefinition as a remedy.We say that ψ satisfies a large magnitude sector constraint, denoted ψ ∼ sect ( q l , q u ) , ifthere exist constants L, M > such that ( ψ ( x ) − q l x ) · ( ψ ( x ) − q u x ) ≥ for all | x | > M ,while | ψ ( x ) | ≤ L | x | for | x | ≤ M . A strict large magnitude sector is defined analogously.This is indeed what happens for ψ ( · ) here, because from Lemma 4 it follows that anychoice q l < − q s < q u will give such a large magnitude sector. Proposition 4.
Suppose ψ satisfies a large magnitude sector constraint ψ ∼ sect ( q l , q u ) with constants M, L . Let c = ( q u + q l ) / , r = ( q u − q r ) / , and suppose the controller K has been tuned such that the loop ( (cid:101) P , K ) is H ∞ -stable and satisfies (cid:107) (cid:101) T ze ( K ) (cid:107) pk _ gn < r − for the peak-to-peak norm. Then for every input w ∈ L ∞ [0 , ∞ ) the non-linear closed loopstate trajectory x cl ( t ) is in L ∞ ([0 , ∞ ) , H ) . Proof:
1) As before let χ = ψ − cI be centered, then | χ ( x ) | ≤ r | x | for all | x | > M , while | χ ( x ) | ≤ ( L + c ) | x | for | x | ≤ M . We show that this implies | χ ( w ) | ∞ ≤ r | w | ∞ + k for someconstant k > and all w ∈ L ∞ [0 , ∞ ) in the time domain. Indeed, sup t> | χ ( w ( t )) | ≤ sup | w ( t ) | >M | χ ( w ( t )) | + sup | w ( t ) |≤ M | χ ( w ( t )) |≤ sup | w ( t ) | >M r | w ( t ) | + sup | w ( t ) |≤ M ( L + c ) | w ( t ) |≤ r | w | ∞ + ( L + c ) M =: r | w | ∞ + k. Note that the same also holds in the truncated version, i.e., | χ ( w ) · [0 ,t ] | ∞ ≤ r | w · [0 ,t ] | ∞ + k for every t > and all w .2) Note that (cid:107) (cid:101) T ze ( K ) (cid:107) pk _ gn < r − means | (cid:101) T ze ( K ) ∗ w | ∞ ≤ ( r − − δ ) | w | ∞ for some small δ > with respect to the time-domain space L ∞ [0 , ∞ ) , and similarly in the truncatedversion.3) But now both (cid:101) T ze ( K ) and the non-linearity χ ( · ) are finite-gain stable in the sensee.g. of [30, Def. 3] with regard to | · | ∞ . Namely | χ ( w ) | ∞ ≤ r | w | ∞ + k and | (cid:101) T ze ( K ) ∗ w | ∞ ≤ ( r − − δ ) | w | ∞ , both fully and in the truncated version. Since r · ( r − − δ ) < , it followsfrom [30, Cor. 1] that the closed loop of Fig. 7 is finite-gain stable in the sense that | z | ∞ ≤ M ( | w | ∞ + | w ψ | ∞ ) + k and | z ψ | ∞ ≤ M ( | w | ∞ + | w ψ | ∞ ) + k for certain M, k > . We derive asbefore that | x t (0 , · ) | ∞ ≤ M ( | w | ∞ + | w χ | ∞ ) + k and | ψ ( x t (0 , · )) | ∞ ≤ M ( | w | ∞ + | w χ | ∞ ) + k for all w ∈ L ∞ , where x t (0 , · ) is with regard to the closed loop.4) Putting Ψ( x cl ( t )) = ψ ( x t (0 , t )) and v ( t ) = Ψ( x cl ( t )) + w ( t ) as before, we can con-sider v ( t ) as a right hand side in the non-homogeneous Cauchy problem ˙ x cl = A cl x cl + v .Accounting for non-zero initial data needs w χ ( t ) = C cl e A cl t x . Since w ∈ L ∞ , we have | v · [0 ,t ] | ∞ ≤ ( M +1)( | w | ∞ + | w χ | ∞ )+ k =: k (cid:48) for all t , and since v is square integrable up totime t , i.e., v · [0 ,t ] ∈ L [0 , t ] , the solution x cl exists on [0 , t ] and is bounded independentlyof t by a constant depending only on k (cid:48) and the decay rate ω ( A cl ) of A cl . This gives x cl ∈ L ∞ as desired, and the solution exists at all times t > . (cid:3) Remark 3.
It is clear that the impact of this result hinges on computing K for a suffi-ciently large sector where the constant k is as small as possible, as that controls how farthe trajectory x cl ( t ) may remove herself from the steady state .6.4. Overshoot.
It has been suggested in the literature that slip-stick is avoided assoon as the non-linear system is globally stabilized. This is obviously misleading, as
ONTROL OF AN ANTI-STABLE WAVE 23 any sufficiently strong disturbance will cause the trajectory x t to attain the value − x t ,however stable the loop. Stability would then only make the difference that the trajectory,after being stuck, returns to steady state when the effect of the disturbance ceases, whilean unstable design might remain stuck. Since the non-linearity ψ ( · ) is concave in theneighborhood of , the term qx t + ψ ( x t ) = ( q + px t ) x t + o( x t ) < qx t is slightly belowthe linearized term qx t , so that a linear controller may overestimate its effect. This maycause overshoot in the response to a disturbance, thereby increasing the risk of slip-stick.That suggests optimizing the closed loop against overshoot in the channel w → y , whichwe realize by simply minimizing the (unweighted) H ∞ -norm of T y w ( K ) . Reduction ofpeak-gain over frequency is known to be a suitable approach for systems with dominantsecond-order characteristics and performs equally well in the present case. In combinationwith the large magnitude sector this leads now to the programminimize (cid:107) T y w ( K ) (cid:107) ∞ subject to (cid:107) (cid:101) T ze ( K ) (cid:107) pk _ gn ≤ /r (cid:107) W u T uw ( K ) (cid:107) ∞ ≤ K ∈ K (29)where T y w ( K ) is the closed loop transfer w → y obtained from plant P , (cid:101) T ze ( K ) refers tothe transfer e → z in plant (cid:101) P , and the channel w → z u in plant P is a safeguard againstunrealistic control actions. This leads to satisfactory results in the blue case, even thoughthe stability certificate is weaker in the sense that the non-linear closed loop trajectory x cl ( t ) is only guaranteed locally exponentially stable and globally bounded. Remark 4.
The peak-gain or peak-to-peak norm (cid:107) · (cid:107) pk _ gn is the time domain L ∞ -operator norm, which for SISO systems is equal to the time-domain L -norm of theimpulse response, or the total variation of the step response [16, Sect. 5.2]. It is harderto compute, let alone to optimize, than the H ∞ -norm, but the bound (cid:107) H (cid:107) ∞ ≤ (cid:107) H (cid:107) pk _ gn is known. Non-smooth analysis of (cid:107) · (cid:107) pk _ gn is beyond the scope of this work and willbe presented elsewhere. In our experiments we use the trapezoidal rule to estimate theintegral of the absolute value of the impulse response of ( (cid:101) P , K ) , and a heuristic to optimizeit. Bounds for (cid:107) · (cid:107) pk _ gn have been discussed e.g. in [10], and a minimization approachvia linear programming is discussed in [19] for the case of full order (unstructured) K .7. Experiments
Gray scenario.
The gray scenario has been addressed with the approach (28), where q l = − . , q u = − . , W u = ss +2e5 . Using Kalman reduction to determine a minimal real-ization, the finite-difference model with N = 50 is used to design a preliminary controller K ∈ K in the class of th -order controllers. The Nyquist test [5, Thm. 1] shows that K already stabilizes the linear infinite dimensional loop exponentially. Moreover, K satisfies the sector constraint (cid:107) (cid:101) T ze ( K ) (cid:107) ∞ = 0 . < r − = 1 / .
64 = 0 . strictly. Afterchoosing a small enough tolerance with (cid:107) (cid:101) T ze ( K ) (cid:107) ∞ + ϑ < r − , we check using [5, Thm.2] that K satisfies even the infinite dimensional sector constraint, so that the non-linearclosed loop ( G nl , K ) is proved globally exponentially stable in the sense of Proposition 3.In a second phase this controller is further optimized with the true infinite dimensionalsystem as described in [5], maintaining the stability and performance certificates alreadyachieved during optimization. Ultimately this leads to the controller K gray ∈ K in (30)which has the same stability certificates, and slightly improved H ∞ -performance. Thiscontroller was then tested in non-linear simulations with spatial discretizations N = 200 .For instance, in Fig. 10 (left) an initial condition θ t (0) < θ t = Ω representing a deviationof from the steady-state was chosen. The controller was switched on at time t = 10 AND DOMINIKUS NOLL and simulated with a square-wave disturbance occurring at t = 15 with magnitude of the steady-state. In the gray scenario linear and non-linear trajectories are almostidentical. That slip-stick may still occur even for this highly stable scenario is seen in Fig.8 (right), but due to stability the trajectory θ t is able to free herself and regain speed. Gray scenario: Rotational speed
Time linear modelnonlinear model
Disturbance
Time -505
Control signal
Time -2002040 linear modelnonlinear model
Gray scenario: Rotational speed
Time -5051015 linear modelnonlinear model
Disturbance
Time -10010
Control signal
Time -1000100 linear modelnonlinear model
Figure 8.
Gray scenario: Occasional slip-stick occurs even with globalstability. Oscillatory disturbance (left). Disturbance at t = 3 , (right). Blue scenario: Rotational speed open loop linear modelnonlinear model disturbance
Blue scenario: Rotational speed open loop disturbance
Figure 9.
Blue scenario: Slip-stick in open loop. K gray = (cid:20) A K B K C K D K (cid:21) = − . − . − . . − . − . .
59 0 0 − . . − . − . . − .
39 0 . − . − . . . . . − . . − . . . − . . . − . − . − K blue = (cid:20) A K B K C K D K (cid:21) = − . − . . . . − . − .
523 0 0 − . . . − . − . − . − . − . − . − . − . . − . − . . . − . . − . − . . − . − − . − (30)7.2. Blue scenario.
The blue scenario is more challenging as the damping parameter λ is between the two critical values λ ( α, q ) < λ < λ ( α, q ) , giving rise to two unstablepoles. Here slip-stick occurs quickly in open loop (Fig. 9). While stabilization of thelinear closed loop is based on the results of section 5, leading to a locally exponentiallystable non-linear closed loop, a global certificate via the sector non-linearity (28) fails due ONTROL OF AN ANTI-STABLE WAVE 25
Gray scenario: Rotational speed
Time linear modelnonlinear model
Disturbance
Time -202
Control signal
Time linear modelnonlinear model
Blue scenario: Rotational speed
Time linear modelnonlinear model
Disturbance
Time -505
Control signal
Time linear modelnonlinear model
Figure 10.
Initial value below steady state, control switched on at t = 10 .Disturbance at t = 15 . Gray left, blue right. Blue scenario: Rotational speed
Time -2002040 linear modelnonlinear model
Disturbance
Time -505
Control signal
Time linear modelnonlinear model
Blue scenario: Rotational speed
Time linear modelnonlinear model
Disturbance
Time -1012
Control signal
Time linear modelnonlinear model
Figure 11.
Blue scenario: slip-stick caused by large disturbances. Stabi-lizing feedback with K blue allows the rotational speed to recover.to the very large primal sector in the blue case. In response, we use the large magnitudesector constraint in tandem with overshoot mitigation. Moreover, a heuristic for thepeak-to-peak norm is used, which leads to the mixed programminimize (cid:107) T y w ( K ) (cid:107) ∞ subject to (cid:107) (cid:101) T ze ( K ) (cid:107) ≤ ρ ( r ) (cid:107) W u T uw ( K ) (cid:107) ∞ ≤ K ∈ K , (31)the parameters now being q l = − , q u = − . and W u ( s ) = ss +2e5 .The idea is to employ the H -norm of the LTI-system in Fig. 7 as an indirect means toreduce (cid:107) (cid:101) T ze ( K ) (cid:107) pk _ gn , which amounts to replacing the L -norm of the impulse responseby its energy. The parameter ρ ( r ) has been estimated using trial and error so that the H constraint ensures satisfaction of the peak-gain constraint in program (29) with parameter r . Starting again with K ∈ K synthesized for a finite-difference model with N = 50 ,we can then certify exponential stability and H ∞ -performance of the infinite-dimensionalloops ( P, K ) and ( (cid:101) P , K ) via [5], and the H -certificate with [4, Lemma 3]. This controller AND DOMINIKUS NOLL is further optimized in the true infinite dimensional system using the method of [5, 4],leading to the final K blue ∈ K in (30). Posterior certification shows that (cid:107) (cid:101) T ze ( K blue ) (cid:107) <ρ ( r ) = 1 . implies (cid:107) (cid:101) T ze ( K blue ) (cid:107) pk _ gn = 0 . < r − = 1 / .
45 = 0 . , whereby thecomplementary large magnitude sector condition is now satisfied in the discretized modelwith N = 200 . Infinite dimensional certification for (cid:107)·(cid:107) pk _ gn is currently not yet available,even though this ought to be established along the lines of [5, Lemma 4, Theorem 3] and[4, Lemma 3]. The controller achieves excellent results in the non-linear simulation. Thisis shown in Fig. 10 (right) where an initial condition generates slip-stick in open loop(yellow area). Triggering control at t = 10 removes slip-stick and additionally providesrejection against strong and sharp disturbances (blue area). Similarly, in Fig. 11 theeffect of switching the controller on is tested on two different disturbances.It should be mentioned that other ways to address the non-linearity ψ have been dis-cussed. In [17] an adaptive controller for a time varying q ( t ) was constructed, while [4]discusses parametric robust control for q ∈ [ q, q ] as well as gain-scheduling of q ( t ) asfurther possibilities. Conclusion
We have presented a novel method to design exponentially stabilizing regulators ofsimple implementable structure for boundary control of a wave equation with non-linearboundary anti-damping. Our results are illustrated in control of torsional vibrations indrilling systems, and two scenarios labeled ’gray’ and ’blue’ are discussed in detail. Weshow that in order to avoid slip-stick it is crucial to optimize H ∞ -performance of theloop. In particular, reducing overshoot by way of H ∞ minimization proved effective forthe more challenging ’blue’ scenario. The ’gray’ scenario had previously been discussedin the literature, and here the substantial improvement of our method over publishedwork is that we can design finite-dimensional exponentially stabilizing controllers, whichin addition show excellent performance. The ’blue’ scenario is new and more challengingdue to inherent instability. We design finite-dimensional controllers which stabilize thewave equation locally exponentially, mitigate the slip-stick effect, and in addition, havea global boundedness certificate, based on the novel concept of a large magnitude sectornon-linearity. References [1] U.J.F. Aarsnes, D. Di Meglio, R.J. Shor. Avoiding stick slip vibrations in drilling through startuptrajectory design.
Journal of Process Control
Control Systems, IEEE , 22(6):2002,37-54.[3] P. Apkarian, M. N. Dao, D. Noll. Parametric robust structured control design,
IEEE Transactioonson Automatic Control
60 (7):2015,1857–1869.[4] P. Apkarian, D. Noll. Boundary control of partial differential equations using frequency domainoptimization techniques.
Systems and Control Letters , to appear.[5] P. Apkarian, D. Noll. Structured H ∞ -control of infinite dimensional systems, Int. J. Robust Nonlin.Control H ∞ synthesis, IEEE Trans. Automat. Control
51 (1) (2006) 71–86(January 2006).[7] P. Apkarian, D. Noll. Nonsmooth optimization for multidisk H ∞ synthesis, European J. of Control
12 (3):2006,229–244.[8] P. Apkarian, D. Noll, L. Ravanbod. Nonsmooth bundle trust-region algorithm with applications torobust stability,
Set-Valued and Variational Analysis
24 (1):2016,115–148.[9] P. Apkarian, D. Noll, L. Ravanbod. Non-smooth optimization for robust control of infinite-dimensional systems,
Set-Valued Var. Anal.
ONTROL OF AN ANTI-STABLE WAVE 27 [10] V. Balakrishnan, S. Boyd. On computing the worst case peak gain of linear systems.
Systems andControl Letters . .[12] M. Barreau, A. Seuret, F. Gouaisbaut, L. Baudouin. Lyapunov stability analysis of a string equationcoupled with an ordinary differential system. IEEE Trans. Automatic Control , 2018.[13] H.I. Basturk. Observer-based boundary control design for the suppression of slip-stick oscillationsin drilling systems with only surface measurements.
J. Dynamic Syst., Measurement, and Control ,139:2017, 104501-1.[14] L. Beji, L. Benchikh. A method of drilling a ground using a robotic arm.
Int. J. Mech. Mechat. Eng.
IEEE Trans. Control. Syst. Tech.
Linear controller design. Limits of performance . Prentice Hall 1991.[17] D. Bresch-Pietri, M. Krstic. Output-feedback adaptive control of a wave PDE with boundary anti-damping,
Automatica
50 (5):2014,1407–1415.[18] N. Challamel. Rock destruction effect on the stability of a drilling structure.
Journal of Sound andVibration
IEEE Transactions on Autom. Contr.
SIAM J. Control Optim. [22] R. F. Curtain, H. Zwart.
An Introduction to Infinite-Dimensional Linear Systems Theory . Vol. 21 ofTexts in Applied Mathematics, Springer-Verlag, 1995 (1995).[23] M.A. Davó, D. Bresch-Pietri, C. Prieur, F. Di Meglio. Stability analysis of a × linear hyperbolicsystem with a sampled-data controller via backstepping method and looped-functionals. IEEE Trans.Autom. Contr. to appear.[24] E. Fridman, S. Mondié, B. Saldivar. Bounds on the response of a drilling pipe model.
Special issueon Time-Delay Systems in: IMA Journal of Mathematical Control & Information , 27:2010, 513-526.[25] Hsiao-Ping Huang, Chung-Tarng Jiang, and Yung-Chen Chao. A new Nyquist test for the stabilityof control systems.
International Journal of Control
IEEE Trans. Autom. Control
44 (1):1999,81–85.[27] E. Navarro-Lopéz, D. Cortés. Avoiding harmful oscillations in a drillstring through dynamical anal-ysis.
Journal of Sound and Vibration,
Springer Proceedings in Mathematics & Statistics
One-Parameter Semigroups for Linear Evolution Equations.
Springer Grad-uate Texts in Mathematics, Springer Verlag, 2000.[30] I.M.Y. Mareels, D.J. Hill. Monotone stability of non-linear feedback systems.
J. Math. Syst., Est.and Control
J. Appl. Anal. [33] C. Roman, D. Bresch-Pietri, C. Prieur, O. Sename. Robustness to in-domain viscous damping of acollocated boundary adaptive feedback law for an anti-damped boundary wave PDE.
IEEE Trans.Autom. Control
IEEE Contr. Syst. Letters
IFAC Proc. Volumes AND DOMINIKUS NOLL [36] D. Salamon, Infinite dimensional linear systems with unbounded control and observation: a func-tional analytic approach, Transactions of the American Mathematical Society
300 (2) (1987) 383–431(1987).[37] D. Salamon, Realization theory in Hilbert space.
Math. Syst. Theory , Toluca, Mexico, 2006.[39] B. Saldivar, S. Mondié, J. Loiseau, V. Rasvan. Suppressing axial-torsional vibrations in drillstrings.
Journal of Control Engineering and Applied Informatics,
SRAIT, 14:2013,3-10.[40] B. Saldivar, S. Mondié, J.C. Ávila Vilchis. The control of drilling vibrations: a coupled PDE-ODEmodeling approach.
Int. J. Appl. Math. Comp. Sci.
Analysis and Control of oil-wellDrilling Vibrations.
A Time-Delay System Approach. Springer Series Advances in Industrial Control,2015.[42] A.F.A. Serrarens, M.J.G. van de Molengraft, J.J. Kok, L. van den Steen. H ∞ control for suppressingstick-slip in oil well drillings. IEEE Control Systems , 18(2):1998,19-30.[43] A. Smyshlyaev, M. Krstic. Boundary control of an anti-stable wave equation with anti-damping onthe uncontrolled boundary,
Systems and Control Letters
Well-Posed Linear Systems , Encyclopedia of Mathematics and its Applications, Cam-bridge University Press, 2005.[45] S. P. Timoshenko.
Vibrations Problems Engineering . Princeton, NJ: D. Van Nostrand Company,1955.[46] R. Tucker, C. Wang. In integrated model for drill-string dynamics.
Journal of Sound and Vibration ,224, 1999,123-165.[47] R. Tucker, C. Wang. Torsional vibration control and cosserat dynamics of a drill-rig assembly.
Mecanica , 33, 2003, 145-161.[48] G. Weiss, R. Rebarber. Dynamic stabilizability of well-posed linear systems. , Miedzyzdroje, Poland, 1:2-9, 1998.[49] G. Zames. On the input-output stability of time-varying non-linear feedback systems. Part I: Con-ditions derived using concepts of loop gain, conicity, and positivity.
IEEE Trans. Autom. Control ,AC-11(2):1966,228-238.[50] G. Zames, P.L. Falb. Stability conditions for systems with monotone and slope-restricted non-linearities.