Capita Selecta

NAVICATIEKUNDE V

college mt 613

CAPITA SELECTA

prof. ir. J . A. Spaans ir.J .H.Wulder Rapport 811-K

INHOUD

-

Herleiden van metingen van scheepsvaste assenstelsel

_naar

lokaal ass ens tel se 1

HELREIDEN VAN METINGEN VAN SCHEEPSASSENSTELSEL NAAR LOKA.AL ASSENSTELSEL

1. Inleiding Prof. ir. J.A. Spaans

Bij een aantal toepassingen op zee is het noodzakelijk am metingen die verricht worden in een scheepsassenstelsel te transformeren naar een lokaal Noord-georiënteerd assensteisel

[fl.

Een Ultra-Short-Baseline (USBL) akoestisch meetsysteem wordt o.a. gebruikt bij het bepalen van de positie van een vaartuig dat met behuip van DP

(dynamisch positioneren) op lokatie wordt gehouden [21.

De sensor van een. USBL systeem bestaat uit een set van drie akoestische transducers die in het viak van het schip zijn aangebracht in langsscheepse en dwarsscheepse richting op afstanden van enkele cm (a cm). Een baken met bekende positie (transponder) op de zeebodem retourneert een akoestische

puls die door n van de scheepstransducers wordt uitgezonden.

In de transducers wordt het faseverschil _4, gemeten in de ontvangst van het signaal in de langsscheepse en dwarsscheepse sensoren. Het golffront wordt viak verondersteld zodat uit de relaties

A$x cos 0

- - .

A/a x 2u lx Cos 0 = . A/a y LiT

de ruimtelij!ke hoeken 0 en 0 voigen waarmee de oriëntatie van het baken vastligt. Celijktijdige meting van de dubbele looptijd bepaalt R zodat de vector van sensor naar baken vastligt in het scheepsvaste assenstelsel (body-stelsel) met Xb-as naar voren Yb-as naar SB en Zb-as naar beneden.

X

Rcos0

b x

= R cos

Zb = R (1-cos2 0 - cos2

e)

= R cos

-1-(1)

Vessel - mounted hydrophone

Projection of R

on to vertical plane

through Y-axis

Acoustic

_{source S}

Projection of

R on to

vertical plane

through X-axis

Fig. 1 Meting van USBL akoestisch systeem in de transponder mode.

Na herleiding van de meting tot het richtingvaste assenstelsel NEV volgt de offset van het vaartuig (AN, AE) t.o.v. het baken. Dit gegeven wordt ver-volgens gebruikt oiu het DP systeem aan te sturen. (Voor duikmoedervaar-tuigen worden drie onafhankelijke plaatsbepalingssystemen voor DP geist!).

Een Short-Baseline (SBL) akoestisch systeem bestaat uit een transponder op de zeebodem of op een volgdoel (bijv. Roy) en drie of vier transducers in het viak van het schip op ruiine afstand van elkaar aangebracht. En van de

-2-transducers zendt pulsen uft waarvan de "Turn around Delay" t. de transponder wordt bepaald in alle transducers.

Bij bekende propagatiesneiheid volgen hieruit: 2 R1

R1 + R2 R1 + R3

} (3) + R4

Hieruit kan de positie van het doel in het scheepsassenstelsel worden vast-geiegd (vier niet-lineaire vergelijkingen met drie onbekenden). Herleiding naar het NEV stelsel geeft offset en diepte van het doel.

Indien de transponder op de zeebodein staat, en is de positie van' de trans-ponder na cailibratie bekend, dan word.t de positie van een referentiepunt op het moedervaartuig berekend uit (3).

Het Artemis systeem is ontwikkeld en wordt gefrabri.ceerd door het

Christiaan Huyghens Laboratorium in Noordwijk. Dit systeem wo.rdt veel toe-gepas:t op DP vaartuigen en op andere vaartuigen die werkzaamheden

ver-richten rond p1atfrms op zee.

Het Artemis syst'eem bestaat uit een FIX-station op. een vast punt (veelal éen platform) en een MOBILE station aan boord van het schip. Beide stations zijn voorzien van een radar-antenne die door een servosysteem op de antenne-as voortdUrend op eikaar ,gericht blijven.

In de MOBILE wordt de "turn around delay" van de EM golven gemeten waaruit de afstand FIX-MOBILE volgt.

Op de FIX wordt d.m.v. een digitale "shaft-encoder" de hoek gemeten tussen een gecalibreerde nuirlchting (bijvoorbeeid Noord) en de richtlng FIX-MOBILE. Deze hoek wordt bitsgewijs overgezonden van dé FIX naàr de MOBILE d.m.v. frekwentie modulatie (FM) op de draaggolf.

Positieberekening vindt plaats op de MOBILE uit de rho-phi gegevens. voor de afstandmeting wordt opgegeven a = 0.75 m en voor de hoekmeting a = 1

boogmin.

In het nieuwe Artemis IlK IV systeem zijn FIX en MOBILE onderling verwissel-baar en wordt ook in de MOBILE een hoek gemeten, en wel de hoek tussen de langsscheepse richting en de stand van de antenne in de. richting van de FIX.

Uit figuur 2 blijkt dat door de metingen e en $ in FIX en MOBILE de koers K van het vaartuig afgeleid kan worden. Deze koers zal aanzienlijk nauwkeuriger zijn d.n die dOor het gyrokompas wordt gemeten.

Een gyrokompas heeft immers een a van tenminste C°.25 terwijl we bier zeker .een. orde beter zijn.

Een nauwkeurige koers is o.a. van belang voor het herleiden van offsets die voligen uit een USBL systeem. Bovendien verhoogt de extra koersmeting de betrouwbaarheid van het totale systeem (Quality assurance!).

-3-Fig. 2 Koers' va'artuigTvoig-t u'it-inetingen 0! en : K = P80° - + e

Voor alie hierboven genoemde metingen is herleiding nodig van de uitgevoer-de meting in het scheepsassenstelsel naar een richtingvast assenstelsel in het horizontale viak.

De meetinst.rumenten die 'hiervoor worden gebruikt zijn de zogenaamde

!Vertjca,1 Reference Uni" (vRU) voor demeting van "pitch" en' "roil" en' bet, gyrokompas voor de meting, van dé kbeshoek.

Een veel gebruikte VRU i's dé Hippy van Dataweil BV in Haarliem. De vertikale richting wordt bepaald door een trage pendulum.

Ala roll-hoek wordtde hoek gemeten tussen de dwarsscheepse Yb-as en het instrumenteie horizontaie viak terwijl ais; pitchhoek d hoek wordt gemeten tusen dë; iangsscheepse Xb-as en bet instrumentele horizontale viak.

Ais het horizontale platform wordt uitgevoerd met een vertikale versnel-iingsmeter die na dubbele integratie de vertikale verplaatsing van ;het vaartuig geeft, dan spreekt men van een "heave-compensator". De output van

dit instrument wordt gebruikt op hydrografische opnemingsvaartuigen om de gemeten waterdiepte te corrigeren voor domp (heave).

Het gyrokompas meet de hoek tussen het Gyro-Noorden en de projectie van de langsscheepse as op bet horizontale viak.

-4-2. Rotaties van Lokaal stelsel naar Body stelsel en omgekeerd

In figuur 3 is een sfeer getekend waarin het scheepsassenstelsel en het

richtingvaste lokale assenstelsel zijn aangegeven. Het richtingvaste

assenstelsel is hier gedefiniëerd als

E naar Oost, N naar Noord en V

vertikaal omhoog, terij1 het scheepsassenstelsel is gedefinieerd als _XB naar SB,

B naar voren en Zb omhoog.

Kiest men de Zb-as omlaag, dan wordt de Xb-as uiteraard naar het voorschip genomen.

In deze paragraaf wordt dus zowel V als Zb "omhoog" genomen.

Fig. 3 Diverse assenstelsels

-5-Bij rotatie-transformaties om de drie opvolgende assen horen de volgende t ransformatiemat rices

1 0 0

R1 (a) = 0 cos a sin a (4)

0 -sin a cos a

cos ci 0 -sin a

R2 (a) = 0 1 0 ₍₅₎

sin a 0 cos a

cos a sin a 0

R3 (a) = -sin a cos ci 0 ₍₆₎

0 0 1

waarbij de relatie tussen de vector x in het oude stelsel en x' in het

nieuwe stelsel luidt:

= R (a) x i = 1, 2, 3 (7)

Een vector x in het richtingvaste ENV stelsel wordt getransfornieerd naar het scheepsstelsel XbYbZb door achtereenvolgens de volgende rotaties uit te voeren, zie figuur 3:

- ENV stelsel naar E'N'V' stelsel door rotatie om V-as over de hoek -K

- E'N'V' stelsel naar xtbyTbztb _{stelsel door rotatie om E'-as over de}

pitchhoek p

- X'bY'bZ'b stelsel naar XbYbZb stelsel door rotatie om de Y'b-as over de hoek r*, zie figuur 3.

Hiermede wordt de totale transformatie:

x' = (r*) R1 (p) R3 (-K) x (8)

Omgekeerd wordt de transforniatie van body- naar lokaal stelsel

x' =R3(K)R1(-p)R2(-r*)x

(9)

waarbij gebruik is gemaakt van R11(a) = R(-x).

De hoek r in (8) en (9) wordt echter niet door de VRU gemeten, we moeten daarom een relatie vinden tussen de gemeten hoek r en de benodigde hoek r*. Hiertoe transformeren we de vector (1,0,0) in het XbYbZb stelsel naar het E'N'V' steistel:

1

R1 (-p) R2 (_r*) 0

0

Omdat de lengte van de vector 1 is, geldt voor de hoek r tussen Xb-as en horizontale viak, waar r negatief gerekend wordt als Xb-as onder de horizon

is:

-sin r = -cos p sin

t*

sin r ofwel: r* = arcsin ( ) ¼cos p cos

t

sin p sin r* -cos p sin r*

-7-(10)

3. Toepassing op de Artemis MR IV

In figuur 4 is N de mobile aan boord en F het snijpunt met de eenheidsbol in de richting van de fix; het hoogteverschil tussen FIX en N is Ah meter. De gemeten peilhoek aan boord t.o.v. het voorschip is ct, terwiji de hoek

rechtsom vanaf het voorschip in het horizontale viak gelijk is aan

_4.

De door de Artemis gemeten afstand M + FIX is de "slant-range" d. Voor correctie van a zijn de hoeken p en r nodig, en het hoogteverschil Ah.

Dit laatste moet worden gechat en eventueel voor geti worden gecorrigeerd.

Een fout in Ah zal overigens alleen op zeer korte afstanden tot de FIX

resulteren in een foutieve correctie.

Fig. 4 Configuratie van meting met Artemis NEC IV

De vector N-FIX wordt bepaald door de beschikbare gegevens: - schuine afstand d

- boo rdhoek a

- elevatiehoek h = arcsin (-) ₍₁₂₎

-8-C

De schijnbare elevatie in het body-stelsel wordt e genoemd. Voor de

benodigde transfortnaties voeren we nog het Artemis-Mobile stelsel XaYaZa in waar de Yaas loopt in de richting M-FIX.

Van het E'N'V' stelsel naar het X aZa stelsel worden achtereenvolgens de volgende rotaties uitgevoerd.

- om de E'-as over de pitchhoek p

- om de Yb-as over de gecorrigeerde roihoek r* - om de Zb-as over de peilhoek - (linksom!) - om de Xaas over de schijnbare elevatie e

N.B. Ult de laatste twee rotaties blijkt, dat de Xaas in de body-horizon ligt.

De transformatie van het XaYaZa stelsel naar E'N'V' geschiedt uiteraard in omgekeerde volgorde met tegengestelde hoeken. Hiermede wordt de eenheids-vector M + FIX in het XaYaZa stelsel na transforinatie naar E'N'V':

(e' R1 (-p) R2 (_r*) R3 (ci) R1 (-e) 1

=

n' (13)

oI

\'

Hieruit is een relatie te bepalen tussen de schijnbare elevatie e en de

ware elevatie h. We schrijven daartoe (13) uit:

&

cos r* sin ci cos e ± sin r* sin e

n' = cos p cos a cos e + sin p sin sin ci cos e - sin p cos r* sin e (14) sin p cos cos e - cos p sin r sin cos e + cos p cos r* sin e

De derde coordinaat is gelijk aan de sinus van de ware elevatiehoek h zo-dat:

sin h = cos e (sin p cos - cos p sin r* sin a) + (cos p cos r*) sin e (15)

Hieruit is e op de bekende manier op te lossen door te stellen:

sin p cos - cos p sin r* sin ci

cospcos r*

tan1

(16)

-9-waarna (15) overgaat in

sin h cos y sin (e + y) =

cos p cos r*

Gebruik makend van (11) wordt (16) dan:

(sin p cos a - sin r

Sifl

a) = arctan

cos p cos r

en volgt e uit:

sin h cos y

e = arcsin ( _{r*) - 1}

De voor pitch, roll en elevatie gecorrigeerde boordpeiling q volgt uit

tan ₌ e' ₍₂₀₎

Uitgeschreven is dat:

= arctan sin a Cos r* + sin r* tan e

cos p cos a + sin p (sin 1* sin a - cos r tan e) (21)

Samenvat tend:

In de Artemis ME( III wordt d in de MOBILE gemeten en e in de FIX. Hoek o wordt doorgegeven aan de MOBILE. De slant-range d kan gecorrigeerd worden voor hoogteverschil d.m.v. r = (d2 - Ah2)+.

De gegevens r en bepalen de offset (AE, N) t.o.v. de FIX.

De Artemis MK IV geeft een extra meting a.

De VRTJ geeft de hoeken p en r, terwiji de ware elevatiehoek h gegeven wordt door (12).

Achtereenvolgens worden de volgende formules toegepast: (11), (18), (19), (21). Hlerinede is de rechterboordhoek ₄ van de FIX vanuit het

MOBILE berekend.

De ware koers K van het vaartuig volgt dan ult K = 180° _- + 0

N.B. Voor "real-time" berekeningen kunnen voor de kleine hoeken p, r, r*,

h, y en e de benaderingen:

-Hiermede wordt: (11) r* -(12) + h =

pCOSa-

rsina

1 - + p2 - + 2 h (1-4 y ) + e -r 2 1-4 p 1 -

4 p2

-(21) ₊ = arctan (1 2) sin a + r* e (1 -

12)

cos a + p (r* sin a - e (1 - + r*

2)

)

De verschilhoek 6 wordt gedefinieerd als het verschil van de gemeten hoek a en de gecorrigeerde hoek :

In fig. 5 is 6 uitgezet als functie van de gemeten hoek a bij een pitchhoek van -1° en elevatiehoeken h van 0°, 1° en 2°. De figuren 6 en 7 spreken voor zich.

In figuur 8 is getoond, dat de correctie voor (pitch = 12°, roll ₌ 0°) en (pitch = roll = 8.5°) verschillend is, hoewel in beide gevallen de Zb-as 12° uit de vertikaal staat. De figuren 9, 10, 11 en 12 spreken voor zich. De figuren 5

_t/m

12 komen uit ref. [3].

(22) sin x = x

PITCH -i ; ROLL S ; H006TE FIX 8 t/ +2 (deg]

Fig. 5

PITCH S ; ROLL -1 ; HOO6TE FIX S t/R +2 (deg]

Fig.

6

12

-848 835 838 825 828 815 818 885 .888 .885 .818 815 828 .825 .538 .835 .845 848 835 838 825 828 815 518 .885 888 .885 .810 .815 .828 .825 038 .535 .848

rt.

lIlihtIL

I ₁₊₈ -8 8

I,

4 9 : 2 I

N

+8,

/

+1__

\\

S -8 8 : 5 8 0

/

\

PITCH -1 ; ROLL -1 ; H006TE FIX 8 tIN +2 (deg]

Pig. 7

Inclinatie in beide qeuallen 12 deg.

8

13

-*1.5 I +1.8 1.5 +2 : :

wr_._.

p-__

_j'

- -: - +

-- t1Tk

TiURW

}ifrL----

'-\

1 e

to

71

7/

_//\N

/1

45 35 ØØ ' 225 2 $ R01L8 315

\\J/

$65 068 855 850 845 040 035 .838 .825 .028 015 818 .885 800 005 .010 .815 .028 .025 .838 .835 .040 .045 .858 .855 .068 865 +8.5 +8.8 8.5

Fir.

g-pig. 10

ROLL +5/

/,4

If

'

7

.. +885

-

:9

25 278

₃₁ 0.98 .'..

\N:

_\

'-....-

/;

PITCH +1 ; ROLL 0

t/i

+5 ; H006TE FIX +1 (deg]

-3 //

/

4/

-5 -

14

-8.00 8.05 8.18

I-

8.15 0.29 I8.25 .0.25 -'-'-'---t9.29 +8. 10 8.10 0.15 0.28

-0.25

0

---..:

+8.20 +0.15 +8.10

/"\\

' \\ +9.05 U, II I!

It

I *

-

_U

/1

315 / /

,f

I ! I

I,

I I PITCH +1 Fig. 12

ROLL 8 tin 45 ; HOOGTE FIX +1 (deq] (elevatie hoogte)

/

RO-5

ROLL +5/'\ / / +4 \

.f

+3

'\

1' . .2 lee 15 -2

I PITCH +1 ; ROLL S t/n -5 ; H006TE FIX +1 (deg] (elevatie

8.10 8.15 8.20 8.25 .8.88 -0. 85 8.18 0.15 8.20 -B.25 hoogte) '.0.18 /. -

-4. Toepassing op het positioneringssysteem van de Benigraph

4.1 Inleiding

De Benigraph, zie figuur 13, is een zeebodem-kaarteringssysteem dat in staat is om een zeer nauwkeurig 3D beeld te geven van de zeebodem.

Hiertoe is een gesleepte "vis" ultgerust met een multibeam sonar die in én "sweep" een strook van de zeebodem opneemt met een bundeihoek van 900

gezien vanuit de sensor in de vis, zie figuur 14. De gebruikte frekwentie is 740 kHz, zodat de vis op niet te grote afstand van de zeebodem kan

opereren i.v.m. het geringe bereik van de gebndkte frekwentie.

Een gevoig van het gebruik van deze frekwentie is bovendien dat niet de

"harde bodem" zal worden gedetecteerd, maar de bovenkant van de sliblaag. Gebruik van een lagere frekwentie zou echter tot onacceptabele afmetingen van de transducer leiden i.v.m. de relatie

b

50A

d

waar b de bundeihoek in graden is, A de golflengte en d de breedte van de vierkante transducer.

Er zijn 256 eleinenten op een plaat van 20 x 20 cm aangebracht met elk een bundelbreedte van 0.75 x 0.75 graad.

- -

_{- -}

AL77TUOE

_-TO CAUATE

APR DATA

Fig. 13 De Benigraph configuratie

- FIPOS

THE APR SYSTt*i

ACOUSTIC F POsITIONNQ

-

16 --S .5--S MILTIBEAM AR (28) PTh SEN _-S.

Fig. 14 De Benigraph multi-beam sonar

Voor een nauwkeurige kaartering van de zeebodem is het noodzakelijk de positie en stand van de sensor continu te weten.

Voor de stand-bepaling is de vis uitgerust met een gestabiliseerd platform voor uitlezing van koershoek, pitch en roll. Het platform wordt tevens gebruikt als Traagheidsnavigatiesysteem (INS).

De door de INS bepaalde positievector en vaartvector wordt in een

Kamanfilter geintegreerd met de vispositie zoals die door een akoestisch systeem wordt gemeten.

In het navolgende zullen we het akoestische vispositioneringssysteem

behandelen.

-Het moedervaartuig is uitgerust met een akoestische transducer. De vis is uitgerust met twee hydrofoons die, symmetrisch t.o.v. bet langscheepse viak, op de staart zijn gemonteerd.

Verder is de vis uitgerust met een druksensor die de hoogte van de

waterkolom boven de vis meet.

Het door de scheepstransducer uitgezonden breedbandige signaal wordt in beide hydrofoons ontvangen. Correlatie van de signalen levert het tijdsver-schil van aankomst At op, zie figuur 15.

.i U LU) N .P.\6_ 1.,,_... a' PUSI dsI.C.

Fig. 15 Correlatiemeting van hydrofoonsignalen.

De scheepstransducer bevindt zich nu op een kegelmantel gezien vanuit de sensor met als centrale as de verbindingslijn van de hydrofoons en halve tophoek:

a = arccos d

waar c = propatiesnelheid d = afstand hydrofoons.

Uit de vertrektijd van het akoestisch signaal van bet moederschip en de aankomsttijd bij de hydrofoons volgt voorts de afstand R.

In figuur 16 is het totale Benigraph systeem getekend. c.At

18

BLOCK DIAGRAM OF THE BENIGRAPH SYSTEM

Fig. 16 Totale blokschema van Benigraph

-4.2, Ret akoestisch positioneringssysteem (APS)

Voor de bepaling van de positie vanhet middeipunt M van de hydrofOons t.o.v. het moedervaartuig zijn de voigende drie metingen beschikbaar:

de vertikale afstand van M tot het wateropperviak verminderd met de !'djepgangt' van de scheepstransducer T; de afstand noemen we d,

de. ruimteiij'ke hoek a tussen verbindingsiij;n hydrofoons en richting M + T,

- de afstand r van M naar T.

We voeren twee assenstelsels in:

Het assensteisel X1Y1Z1 met de oorsprong M, de Z1-as vertikaal omhoog,, de Y1-as evenwijdig met de voorliggende koers (heading) van het moeder-schip en de X1-as dus naar SB in het viak door M evenwijdig met het horizontale viak.

Het body-assentesel XbYbZb met H ais oorsprong, Zb-as totnhoogtt en Yb-as in: het symmetrieviak van de via, naar voren.

De drie meetvergeiijikingen zuilen worden uitgeschreven in bet

body-assen-steisel,, waarbij de index b wordt weggelaten. We zoeken dus de piaats

vector (X,,Y,Z) van T in het body-stelsel [4].

De X-coordinaat van T in bet body-steisel is de projectie van de afstand r op de Xb-as, waar de hoek tussen de Xb-as en MT geiijk is aan a:

x = r cos (a) ₍₃₀₎

De X-coOrdinaat is hiermede direct te berekenen. De gemeten afstand r levert op:

Y + Z

= r - X = r sin a

2 .2 2 2 2 2

(31)

De diepte d moet worden getransformeerd naar bet body-steisel. .Dë vector (o,

0,

d)

= wordt voigens (8) getransformeerd naa,r bet body-stelsel. AIs voor hoek K het verschii genomen wordt van de koersafiezingen van moederschlp en vis dan moet in (8) +K ingeiu1d worden i.p.v. -K.

(r*) R1 (p) R3 (K) d (32)

-Ms

de hoek is tussen de gevraagde vector en de vector dan geidt voor het inwendig prOduct:

= iikit

iiii

cos = d2., zodat dc derde vergelijking.voor de

oplossing van (X,Y,Z) luidt:

X Xd + + Z Z1 = d (33)

Dit is de vergelijking van een plat viak ioodrecht op dé normaaivector (Xd,Yd;,Zd)b = (o, o, d)1 door de transducer T.

d2 - X X - ZdZ

Uit (33) volgt Y = d A - BZ (34)

d

Dit substitueren in (31) geeft

Z2 (1 + B2) - 2 ABZ ± A2 - X2 - = C (35)

Hieruit wordt Z o:pgeiost, waarbij van de opiossingen ais onteaiiistlsch afvalt.

De gevonden waarde van Z ingevuide in (34) geeft Y.

Vervoigens wordt de _{gevonden (X,Y,Z)b} teruggetransformeerd naar het 1-steisel:

= R (-K) R1 (-p) R2 (_r*) (36)

zodat de positie van het punt T t.o.v. M in het 1-steisel vastligt. De

positie van H in het naar T getransleerde i-steisei is dan

-[11 Spaans, J.A.

Note on coordinate transformations

Internal note Intersite Surveys G 473/JAS/ISC Haariem 1986

Hilne, P.H.

Underwater acoustic positioning systems E. & F.N. Spon Ltd London. 1983

Wardenier, R.

Het Artemis .plaatsbepalingssysteem

Afstudeerscriptle HZS A'dam afd. Hydrografie 198W

ter Steege, J.W.

De positionering van de Benigraph

afstudeerscriptie HZS A'dam afd:. Hydrograf.ie 1988

-DYNAMISCH POSITIONEREN

Johan Wulder

In dit college wordt een overzicht

_{gegeven van de factoren die}

een rol spelen bij het ontwerpen van een dynamisch positionerings

systeem. Hierbij zal gebruik worden gemaakt van de kennis welke

opgedaan is in de colleges mt612 en w61.

1 inleiding

In het college mt612 is een methode

_{gegeven om met be.hulp van het}

kalman filter de positie van

een vaartuig te schatten. De stap

die hierop volgt is

_{er voor te zorgen dat men zich niet op een}

gescha.tte maar dat men zich automatisch

_{op de gewenste positie}

bevindt. Deze gewenste positie kan é.en stationair

_{punt zijn of}

een geplande baan. in dit college wordt van

een stationair punt

uitgegaan.

Dit

geval doet zich voor bij

o.a.

_boorschepen

_en

duiker support schepen. Ret automatisch volgen

van een geplande

b:aan is van belang voor bijvoorbeeld seismisch

_onderzoe.k.

Men kan een vaartuig op twee manieren

op

een gewenste positie

houden. Ten eerste door een vorm

_{van verankering, ten tweede door}

g:ebruik te maken van de voorst.uwing

_{van het vaartuig. Dit laatste}

heet 'dynamisch positioneren'

.

Ook combinaties van beide methoden

komen voor.

2 Het DP-systeem

Uit de inleiding volgt dat bet doel

_{van een DP-systeem is het op}

zijn plaats houden van een vaartuig met behulp van zijn

voortstu-wing. Het isduidel.ijk dat dit brandstof (=geld)

_kost.

Verder wil men ec'hter ook dat bet schip niet

_{te veel van zijn}

gewenste positie afwijkt, bijv.

_{in verband met .de riser van}

_een

boorschip.

Met deze twee

tegenstrijdige

eisen

zal

rekening

gehouden moeten worden. Een definitie van dynamlsch positioneren

kan zijn: Een vaartuig, voor bepaalde

_{ontwerp kondities, binnen}

bepaalde toleranties op een gewenste positie houden met minimale

brandstof kosten.

Om bovenstaande te bereiken moet

_{een regelaar ontworpen worden}

welke continu een afweging uitvoert

_{tussen de positie afwijking}

en de stuwkracht. In figuur 1

_{is de plaats van deze regelaar in}

Xg ,Yg ,

dat cie

+ ,

ge:en

tijd

Normaal

is

maken van een

Deze aanname

regelaar.

Verder heeft

RE GE LAAR

Fstuw

SCHIP

w(t)

-X,Y,!/

toestand van hat schip: In dit college wordt

_{er van uit gegaan}

positie en de sneiheid direct gemeten worden

_{en dat er}

-

verloren

gaat met

het

verwerken van het

signaa.1.

dit niet bet geval en zal

men b.v.

gebruik moeten

kalmanfilter om de actuele positie

t:e voorspellen.

heeft verder geen invloed

_{op het ontwerp van de}

men wel

de sy.steem verg.eilijkingen van he-t.-.schip

met: Xg,Yg,.g: gewenste positie en koers

x,Y,.p

_{werkelijke positie en koers}

w(t) s toringen

Fstuw

_{gevraagde stuwkracht (vector)}

figuur 1

Regelkring Dynamisch Positionerings systeem

Indien men deze regelkring bekijkt bl.ijkt

dat

er naast de

te

ontwerpen regelaar nog drie

'parameters'

een rol

spelen. Deze

zijn:

de storingen

toestand van bet schip (positie en koers)

voortstuwing

Wanneer men een regelaar wilt ontwerpen zal

_{men van deze factoren}

enige kennis nodig hebben. Dus:

-

storingen:

Deze worden voornamelijk door de omgeving

veroor-zaakt,

namelijk:

wind,

stroom

en golven. De

grootte

van de

krachten

i.s

afhankelijk

_{van het}

_gebied

_waarin

_het

vaartuig

operationeel

is. Bij het

ontwerp

wordt

dit

meestal

door de

opdrachtgever bepaald.

_{Deze geeft de maximale windkracht en/of}

g'olfhoogte

op waaronder het

systeeni moat functioneren.

_{Van de}

dri

_{genoemde storingen kan men het vo.lgende}

_{zeggen. De wind en}

de stroom geven elk een kracht

_{en moment die ongelijk aan nul}

zijn.

_{De golfkrachten veroorzaken drie}

_{soorten bewegingen.}

_Ten

eerste hoogfrequente bewegingen welke dezeifde frequente

hebben

a.ls de goiven, met een gemiddeide van nul. Daze zijn voor een

DP-systeem

niet

interess.ant

omdat

ze een zeer

geringe

positie

afwijking veroorzaken die bovendien gemiddeid

_{nul is. Ten tweede}

de laag frequent.e bewegingei-i

_{welke veroorzaakt worden door niet}

iineaire effecten in het goifveId. Deze maken

_{bij stationaire DP}

systemen ongeveer 70 procent van de stoorkrachten

_{uit. Ook deze}

veroorzaken een gemiddelde beweging

_{van nul. Tot slot veroorza.ke.n}

de golven een gemiddelde driftkracht. Hiermee

_{zijn de}

beiangrijk-ste stoorkrachten beschreven.

met M

:

massa vaartuig

Fx: stuwkraàht langs x-as

B :

demping term

Fy: stuw.kracht lan.gs y-as

M :

moment veroorzaakt door stuwkracht

Mw: moment veroorzaakt door storing

Fwx en Fwy

:

stoorkrachten langs respectievelijk de x en y-as

Indien men deze vergelijkingen in de toestand

_{vorm schrijft, dan}

voigt voor de langsscheepse vergelijking:

met: Xm

xl x2

Voor

de

andere

verge hf kingen.

Men

kan

deze

voigens:

Xm

[m2I

k 1

m/b(1et/m)

0

et/m

de gemeten positie en sneiheid.

x

x

variabele.n y en b

gelden

dezelfde

_toestand

vergehij kingen

1Ix11

I I +

v(t)

±1 Lx2Jk

F/b+m/b2ebt/mm/b2

L/b(1ebt/m)

de

discrete

vorm

schrijven

Fx-f-(t)

(7)

(8)

nodig. Deze zijn t.o.v. een aardvast assenkruis:

Mit + B

Fx + Fwx

(1)

MI + Bt

=

Fy + Fwy

(2)

IZZjb+B,.=M+Mw

(3) + w(t) + w(t) (nit6l2 paragraaf

(4)

(5) 10.3) k+1

+ f(t)c(t)dt

met eFt (6)

(zie blz 65 mt612)

ofwel: 0 Fwx/M 1 0 xl +

x2

1 o 1

To

L'

o -B/M

x2

R

ol

I

II

I

+v(t)

Xm2

[

1

Li

ofwel: x(t) = F (t) + C u(t)

xm

H (t) + v(t)

of anders geschreven:

.m =

Hxk

-

voortstuwing:

In dit colleges wordt er van uit

_{gegaan dat de}

gevraagde stuwkracht gelijk is

aan de geleverde stuwk.racht.

In

werkelijkheid zal een apart programma de thrusters

_{aanwijzen die}

de

stuwkracht

moeten

leveren. Dit

gebeurt

omdat

bij somm:ige

combinaties

van

thrusters een

aanzienlijke

rendement

daling

optreedt

door

interactie van

de versch:illende zogstromen. Dc

combinaties van thrusters die gebruik;t worden is

zeer ve!rschil

lend,, bijvoorbeeld:

Voith Scheider

normale schroef en boegschroef

langs een verticale as roteerbare schro:even

etc.

3 Ontwerpen van de regelaar

Zoals al opgemerkt moot dé regelaar

_{twee teg:enstrijdige eisen met}

elkaar verenigen, namelijk minimaal bran'dstofverbruik en minimale

af:wijking.

Er bestaan verschillende s:oorten rege.laars waaruit men een keuze

moet make.n,

namelijk de

P, P1,

PD en de PID regelaar.

Elk van

deze heeft zo zijn eigen voor- en nadelen. De P regelaar werkt

alleen storingen met een gemiddelde van nul we.g. Wanneer dit niet

het

geval is zal

er

_{een statische afwijking ontstaan.}

_{De P1}

regel.aar heeft van het voorgaande

_{geen last echter men kan een}

onstabiel

p.roces krijgen. De PD

regelaar

geeft een snelle

responsie en kost dus veel energie. De PID

_{regelaar combineert de}

eigenschappen van bovenstaande re.gelaars.

Voor de

regelaar van

Oris

DP systeem wordt uitge.gaan

_{van een P}

regelaar

(

welk p:robleem verwacht je

?).

_{De afleiding voor de}

andere type reg:elaars gaat op de zelfde wijze.

De eerste stap bij het ontwerp

_{van een regelaar is dat men het}

systeem lineariseert rond het werkpunt. Vervoigens

_{gaat men de}

afwijkingen ten opzichte van dit werkpunt

_{wegregelen. Ons systeem}

is al lineair we moeten riu aileen

_{nog een werkpunt kiezen en de}

afwijkingen

toy dit punt

beschrijven.

Dit

doen

we door de

oorsprong

van

het

assenkruis

In de

gewenste

positie

en de

richting van de x-as evenwijdig

_{aan de gewenste koers te kiezen.}

Er tre.den dan enige veranderi.ngen

_{op toy figuur 1.}

_{De ingang van}

het. systeem is nul en de uitgang is niet

_{langer de positie niaar}

de posltie afwijking. In figuur 2

s dit nog eens getekend.

26-rege.laar

figuur 2

regelkring toy werkpunt

De eerder afgeleide formules gelden

_{nog steeds met dit verschil:}

in de toestandsvector staat niet langer de vector X= [X

]t _maar

de afwijking toy deze vector:

We willen flu u(t) en x(t) minimaliseren volgens:

J =

1/2E(4Q1x

_{+ Q2Mk )}

is minimaal

(9)

onder de voorwaarden

Zk+1

+ruk

(4)

xmk

=

H2ck

(5)

Dit gaat volgens methode met

_{de Lagrange multiplicatoren}

_(zie

mt612 blz3O).

De te minimaliseren functie wordt hierdoor:

k=0'

1/24Q1k + 1/2uQ2uk +

k+1(k+1k+Tk)]

(10)

De vergelijking (5) wordt straks in rekening

_gebracht.

De functie J moet flu geminimaliseerd warden naar uk,xk en

Dit levert de volgeride drie vergelijk-ingen

_op:

='uQ2+A1r=0

6J k+I

-k+1+k+rk

0 ₍₁₂₎ 6J

.Qi -

(13)

Indien

en _o

_{bekend zouden zijn kan voor elk tijdstip}

k

berekend worden.

_{is echter onbekend. Wat we wel weten is dat}

u=O, omdat deze geen invloed

_{meer heeft op}

_x

_{en men de}

_eis

heeft gesteld voor minimaal energie

_{verbruik. Uit (11) volgt dan}

waarna met (13) volgt:

= Qi (14)

5stuw

) Schip

Hiermee

is bet

probleem gedeflnieerd.

_{De uitwerking van deze}

vergelijkingen is in bijiage 1

_{gegeven. Het resultaat iuidt:}

11k -KkXk (15)

met

_{Kk (Q2rtsk+lr)l rtsk+l}

₍₁₆₎

met _Sk=

tMkl

_{+ Qi (17)}

en

_{Mk+l= sk+1.skir(Q2+rtsk+lr)lrtsk+lg}

_(.18)

met de randvoorwaarde:

_{s= Qi}

₍₁₉₎

Wat houdt dit in

?

Men heeft als randvoorwaarde

_{S=Qi. Men zal}

dus op het tijdstip n uloeten beginnen met bet berekenen van K en

vervoigens

_{terugrekenen tot het tijdstip kO. Ondertussen}

_moet

men de versterkingsmatrix Kk voor k=0,...,n opslaan om later te

kunnen gebruiken.

Welke prob.lemen geeft dit in de praktijk ?

Ten eerste heeft men

een groot geheugen ruimte no:dig om alie versterkingsmatrices

_op

te kunnen slaan. Ten tweede moet men het aantal tijdstappen

_weten

welke men wilt regelen. In de praktijk is dit

_{niet bekend.}

Er

blijkt

gelukkig

een andere

mogelijkheid

te zijn, welke

suboptimaai

s.

_{Uit de theorie blijkt dat de versterkingsmatrix}

tijdsafhankelijk is. In de praktijk is K ec'hter constant over het

groot.ste

deel van

de tijd

en aileen

de

iaatste

regelstappen

kenmerken zich door een afwijking.(waarom ?)

Omdat we bij dynamisch positioneren aitijd

_{ee.n t1oneindige" tijd}

willen regelen zijn we alleen geinteresseerd

_{in het constante}

deel van de versterkingsmatrix. Ofwei:

_{we gebruiken K0 over de}

gehele regelperiode. De voordelen zijn duidelijk.

_{Bij het ontwerp}

van

de

_{regelaar doorloop je eenmaal bovengenoemde algorithm.e,}

waarbij je ailee.n K0 hoeft te onthouden.

Het

laaste

onderdeel

welke we

nog moeten definieren zijn

de

weegmatrices Ql en Q2 Deze worden tijdens het

ontwerpen van de

regelaar bepaald en zijn afhankelijk

_{van de verhouding tussen de}

stoorkrachten en de maximale positie afwijking.

_{In de volgende}

paragraaf wordt hierop teruggekoinen.

_{De stru'ctuur van}

_{Ql is wel}

bepaald.

D:eze is

_{bepaald door de nog niet gebruikte}

randvoor-waarde

(5). We

hebben met metingen

te

maken terwijl

in de

afleiding van de complete state vector is uitgegaan

Pus

we hebben gebruikt

: ipv

mtQ1 x

Echter met

_'m

₌ H& volgt:

waar uit voigt

HtQ1H

-

-26'-Q IxmtQl.mI_I:XtHtQlHxI

JLZ

10 1

Lc1

0 0 c2

met c: weegfactor

Nu is de regelkring geheel gedefinieerd

_{en kan het schip op zijn}

4

Simulatie

In

de

vorige

_paragraaf

_is

_de

_{versterkingsrnatrix afgeleid.}

_In

principe zou je hiermee direct

_{een schi,p kunnen regelen. In de}

praktijk is dit te riskant en test

men een regeling door middel

van

een

simulatie.

De

tweede

reden

_om

gebruik

te

_{maken van}

simulatie

technieken zijn

de

weegniatrices Qi en Q2 Met deze

factoren kunnen

de

maximale

_{positieafwijking}

_{en het vermogen}

gewogen worden. Tijdens de siniulatie zal Inoeten blijke.n of deze

interdaad

niet

buiten

de

maximale

waarden

vallen.

Met

het

criteriurn J wordt niets over de maximale waarde gezegd! Stel dat

je de maximale verstoringen kent

_{en je rnerkt tijdens het}

sirnule-ren dat je een grotere positie afwijking krijgt dan toegestaan,

dan kan je door niiddel van de weegmatrices deze afwijking kleiner

maken. Het benodigde vermogen neemt dan echter wel

toe

1

Indien

men nog moet beslissen welk verrnogen er in het vaartuig

g:einsta1-leerd moet

worden kan

men

flu

uit

de

simulatie het

maxirnaal

benodigde

vermogen

_(-

het

minimaal

te

installeren

ver.rnogen)

bepalen. Hieruit volgt dat, naast het

_{testen van de regelkring,}

een simulatie ook kan dienen als een stuk gereedschap

_{bij het}

ontwerp van een DP-systeem.

Als men van een simulatie

_{programma gebruik maakt zal men zich}

goed moeten realiseren dat

men met een model van de werkelijkheid

bezig

is.

Er

zullen dus

in

de

werkeiijkheid altijd effecten

optreden welke niet gemodelleerd zijn.

_{Bij het gebruik van een}

bestaand model

zal

men

altijd

moeten proberen

_de

vereenvou-digingen in het model te achterhalen.

Voor dit college is een (eenvoudig) simulatie

_{model opgesteld van}

een

DP-systeern.

Hierbij

is

van

de

_{vergelijkingen}

_1,2

_en

₃

uitgegaan, waarbij:

M

=10000

kg

=2000

kg/s

B

=4000

kg/s

Izz-50000.0 kg/rn2

B

=5000

kgm2/s

De versterkings matrix K is met (l5)-(18) berekend, resultaat:

K= [

1 _9.951

5= [

1 4.98]

K -

[ 1 98.9]

w(t)

K l/M B

----H'

I-Figuur 3

blokschema voor MX+BX

_xwx

Vervolgens is met de simulatie taal "TUTSIM" het

_systeem

gesimu-leerd. Ais voorbeeld is vergelijking (1) als blokschema

_getekend

in

figuur

3.

_{Hieruit blijkt dat interacties tussen de diverse}

bewegingen verwaarloosd zijn. Het complete schema

_{is in bijiage 2}

gegeven. De stoorkrachten worden nagebootst door een sommatie

van

vijf sinusgolven. In de figuren 4

_{en verder zijn enige resultaten}

gegeven voor:

-

x,y als functie van de tijd

-

als functie van de tijd

-

x,y met stroom als functie van

-

de beweging in het x,y,vlak

Andere combinatie.s zijn natuurlijk ook mogelijk.

_{Het zal}

duide-lijk

zijn dat

simulatie

een krachtig

gereedschap

is bij het

ontwerpen van een regelaar.

LITERATUUR

Lewis F.L.

Optimal Control

1986, John Wiley and Sons

ISBN 0-471-81240-4

Morgen M.J.

Dynamic positioning of offshore vessels

PPC Books, Tulsa, OklahOma,l978

Olsder G.J.

Filter en identjfjcatjetheorje

college a42

februari 1986, TU Delft, Wiskunde en Informatica

30

Ilod(21.

[j

] .:

Date

9

/

:[5

/

:19fl7 11fl1E

6

17

Ti

nu

nq

o

i000000

, DELTA

200,0000

,RANIE

ic:'t1

id

otI31ock

F:L

For

L.1. oc:kNo .,

Fl

ot-M

INi

mum , Plc:rL-MAX.i mum Comrnrut. 0 I-10r2 ,

0.

0000

.

200

0000

ii

m.

vi

.. , --100,,

0000

, 100.. 000(1 Y2 8 .1-100..

0000

., 100..

0000

1.3 Y3

,--

[00..

opoo

,

i.00.

ooc:o

Y4:

:19 ,

--:[0..(00()

,

--2.

c)000 Ut' £.Lu 1U Co uQfl .uv'UvvI

4I

1)Q' nwo IIfI VVW L I I I I I I I F1GUUR 4

plot

x,y,

-tijd

_3,-

I I I I I I. I

I

NCJdE1 1 JE

dp

Datc

9 /

:15

/

1987

Tiinc

8 7

1' 1 mi.

0. 5000000

,

P1. c:jtEU'cjc:ks and 8:ai cs

F ci mat

Lti:ocI::I\Io,

_Plot-MINI

mttiI1,

J-!or

0

, 0

Yl

,sc)c) c)c)(:)c)

8

, -500 - 0000

1 9

C00

(flW 'IcIvIUt,w.

221t33

fM

£VU, U1!VU n OQflhl

-5@ 12a

S 000E+03.

, RANGE P1 o1:--MXirnum Comment 5. c:x:)oEc):3

Ti

SOC). 0000 :,(_)() ()(..)(_)(.)

--2 0000

--:... t)c)c)(:1 x1-t.p1:ot stPoaH I I I I J J I

Jij VViV/i'JVvvvr

WMM!WVVV\/V\MWMWMWMAAM

L L H)flfl UtUU!J

7-FIGUUR 5

plot x,y-tijd met stroom

5;1f3

i::c?J.

I:I i.

Datp

9 /

15 /

:L9a7 .

1img

b_4

Ti ml n

0.

5C:oc:Ku:J , DELTi J. c:x:)cJE+c:)3

,RiNGE

P:LotB:I.ocl::s. and

3ca:Le

For mat

E']. ockNo

RI. ot-MINi mum.,

I-]. ot-MAXi rnum Comment

Hor

,

0. :)o0o

:1

000E+03

Ti me

::

,, -6ico 000)

2O.O 00.00

Y2 :38 ,

--10 000E0:3

,

:;o,. 000E+03

-jo.. 0000

.,

-2. 0000

Y4 :19 ,

-10 0000

-2.

0000

_1fl flIlW) JiIvI.v1IUv

28Lt3

-368 89 -44 -6gL 1 (U% A II ZWA I 1 . 1 I I I I I

L3

iie

ii3

-23-FIGUUR 6

plot Fwx,x-tijd

f'

'[F

Dat

9 /

% /

1987

Time

7

7

I1 ifi L flr4 ( ) L,( )( )( )(it IC I DEL I

., J IC It F+ r

hNrir

FlotBlo:ks, arid Sca1es

Format

6 481.g@ _*I 171]fifl 1W, !JUVCI I I .1 I I I I

-FIGUUR 7

plot x-y

I

xi plot

e.t

stMM

1u

I.,

Fl cit-MINI mum

.,,

P1 otHiXi mumi

Comment

Hor

- 100. 0000

_,

100 - 0000

8

--i. 00., 0000

,., 100. 0000

Y2 Y4

8

8

., 1 9 _,

-'1oo., 0000

_,

- i. 00. 0000'

_,

-1 0. 0000

,

:ioo. 0000

100. 0(100

-2 000c)

OPTIMAL CONTROL 251

9.3 OPTIMAL CONTROL

Optimal control methods are attractive because they handle multi-input systems easily and allow the designer quickly to determine many good candidate values of the K matnx We will develop the time varying optimal control solution first and then reduceit to a.steadystate solution. This amounts to another method of corn-puting the K-matrix in the control law

u=Kx

(6.19)

discussed in Chapter 6 and in the examples of Section 9.2 (EqS. (9.2) and (9.3)].

9.3

Tlme-Varyln9 Optimal Solution

Given a discrete plant

Xk+1 4Xk + r'uk, (9.10)

wewish to pick Uk so that a 'cost function'

+ UrQ2Uft] (9.11)

is minimized. Q and Q2 are symmetric weighting matrices to be selected by the designer who bases his choice on the relative importance of the various states and controls Some weight will almost always be selected for the control (1Q21 * 0)

252 MULTI VARIAB AND OPTIMAL ONTRUL

states would be dnven to zero at a ridiculously fast rate which could saturate the actuatordevice.2 The 0's must also be non-negative.3

Another way of stating the problem gi'ien by Eqs. (9.10) and (9.11) is that we wish to minimize

I

_f

[xQ1x, + UIQSUft] k-a

sUbject to the constraint that

Xk+l+4Zk+FU.O,

k=Q1, ...,N.

(9.10)

Thisis a standard const ed-il nima problem which can be solved using the method of Lagrange multipliers There will be one Lagrange multiplier vector,

A11,for each value of k. The procedure is to rewrite (9.10) and (9.11) as:

NI

I

=

_{xQx +}

+ Ar+z(xk+, + Xk +

rukJ

(9.12)

and find the minimum off with respect to x,Uft,andAk. The subscript on A is

arbitrary conceptually but we let it be k + I because this choice will yield a par ticula1y easy form of the equations later on.

Proceeding with the minimization leadS to:

uQ3,+ AI+IT = 0. control equations, (9.13)

x1

+

+ I'ti = 0,

State equations, (9.10)

and

= xrQ1 - xr +

g(I) = 0.

adjoint equations. (9.14)

The last set of the. equations, the adjoint equations, may be written as

Ak = (9.14)

Combining (9.10), (9.13). and (9.14) gives uS a set of coupled difference equations

defining the optimal solution of Xk and A and hence Uk provided the initial (or

final) conditions are known The initial conditions on x must be given however

usually A would not be known, and we are led to the endpoint to establish a final condition. From Eq. (9.11) we see that UN will be zero for the minimum 5 because

'11 the sáthlin rate. T. is long. hówevcr.,m control which moves the state aiongas rapidly as possible

may be feasible Such controls are called dead beat since they beat the state to a dead stop in at most n steps. They correspond to placement of all poles at Z = 0.

Matrix equivalent of a nonnegative number it ensures that2TQ,x and UTQU are nonnegative for afl

possible a and a.

(9.11)

ar U UU,Ia..-wnUn5p J

11N has no effect on XN[see Eq. (9.10)]. Thus Eq. (9.13) suggests that AN+a = 0,

and Eq. (9.14) thus shows that a suitable conditiOn is

AN = QIXvf (9.15)

The solution to the optimal, control problem is now completely specified. It

consists of the two difference equations (9 10) and (9 14) with u given by (9 13) the final condition on A given by Eq (9 15) and the initial condition onx0assumed

given in the problem statement. The solution to this two-point boUndary-value problem is not so easy.

One method called the sweep method by Bryson and Ho (1969) is to assume

= (9.16)

then (9.13) becomes

Qsua =

=

P$k+i(4'xa +

Solving for u, we obtain

Uk = + 1'Sft+ifltFTSk+lVxk

= R1'Sk+id'xa.

(9.17)

In (9.17) we have defined R = Q + rs+r' for convenience. If we now

substi-tute (9 16) into (9 14) forA5and A5 we eliminate A Then we substitute (9 17)

into (9.10) to eliminate Xk,.i as follows. From (9.14)

Ak = 4STA _{+ Qixi}

and substituting (9.16), we have

SkXk = TSsXa + Q1x1.

Now we use (9.10) for x54.1,

Skxk = 4TS5+1(4'x + ru5) + Q1x5.. (9.18)

Next we. use (9.17) for u in (918),

S5Xft = 43TS

(tlx5 - fRtI'TSk+ixk) +

and collect all terms on one side:

S5 .S+1 +

Sk+tfR S5 - Q1]x5 = 0. (9.19)

Since Eq. (9.19) must hold for any x5 whatever, the coeffident matrix must be identically zerO, from which follows a backward equation in Sk,

= T[Sk _{- Sft1FRtI'1Sk+J} + _Qi, (9.20)

which is often rewritten as

254 MULT1VARIABLE AND OPTIMAL CONTROL

where

= Sk+i- Sk+iflQI + PSk+ifltI'TSk+s. (9.22).

Note that the matrix to be inverted (R) has the same dimension as the number of

controls which is usually less than the number of states

The boundary condition on the recursion relationship for S1+1 is obtained from (9.15) and (9.16); thus

SN = Qi. (9.23)

and we see that Eqs. (9.21) and (9.22) mustbe solved backward with the iflitial

condition (9.23). To solve for u, We use (9.17) to obtain

Uk= -Kkx,

(9.24)

For any given initial conditiôh for x, to apply the control, we use the stored gains

and

= 4iXk + Fuk, (9.10) where

= -KkXk. (9.24)

Notethat the optimal gain.Kk, changes at each time step but can be precornputed

and stored for later use as long as the length,N. of the problem is known; for

ex-ample. no knowledge of the initial state x is required for computation of the con-trol gain Ka.

As an example of the time-varying nature of the control gains. Eqs (9.21)

through (9 25) were solved for the satellite attitude control example described in

Appendix A.1. The system transfer function is I/s2.

The state weighting matrix was arbitrarily chosen to be

11 0

Qi-o

_o'

which means that the ang1 state was weighted but not the angular velocity. The control weighting matrix is a scalar in this case because there is a single control input; three values were selected for evaluation:

Q2=0.01,0.1,and1.0. (9.26b)

The problem length for purposes of definingwas chosen to be SI steps. vhich, with the sample period of T = 0.1 sec means that the total time was 5.1 Sec.

Figure 9 3 contains the resulting gain time histories plotted by the computer We see from the figure that the problem length only affects the values of K near

the end and in fact, the first portions of all cases show constant values of the

0.5-00 - I I

0 10 20 30 40 50

Fig. 9.3 Example of control gains vs time

OPTIMAi. CONTROL 255 (9.26a) 9 8 7 6 C C₀ 5 2 4 C 0 U 3 2 0 4.5-4.0 3.5- 3.0-C C 2.5-2 C 2.0-0 U 1.5-1.0 Q,=O.OI Q2 =0.1

-I.0

I I 10 20 - -- 30 = 0.01 - 40 50 =0.1 Q2 1.0 = (Q2 +

and is the desired "optimal" tithe-varying feedback gain. Let us now summarize the entire procedure:

(9.25) 'I I.

LetS,-QlandKN-O.

(9.23) 2. Let

k4-N.

3.

LetMft._Sk.Sk+1'1'SkilTTSk.

(9.22) 4. Let K,

- (Q + F%FISkD.

(9.25)

.IJ

5. Store Kk_I.

I 6. 7.

LetSk-47M+Q1.

Letk4-k-l.

(9.21) 8. Go to step 3.

gains. Thismeans that for problems of long duratiOn (longcomparedtotransients

:in the gains), the optimal controller for much of the time is identicalin structure to

the constant gaincasesdiscussed in Chapter 6 and Section 9.2 except that the

val-ues oftheconstant gainsarebased on a cost function rather than root locations.

For multi Input problems the time-varying gains act exactly like the example above The next section develops a method to compute the constant value of the optimal gains so that they can be used in place of the time varying values thus yielding a much simpler implementation and ye! one that is almost optimum.

Before we leave the time-varying case it is informative to evaluate the op-timal cost fUnction in terms of A and S. If we substitute (9.13) and (9.14) for

g.i' and A+'& in (9.12), we find

f =

[4Q1 + UQUE Ar+IX+I

+ (Al - xlQi)x + ( t)u]

= _{(AIX Al+lxk+l)}

Axo-1AN+jxN+i.

However, from (9.15), AN+t 0, and thus, using (9.16), we find

=I4SGXO.

Thus we see that having computed S; we can immediately evaluate the cost asso-ciated with the control.

9 3 2

Steady-State

Optimal Solution

The optimal solution for the time-varying gain. K. is described by Eqs. (9.13), (9.14), and (9.10). We can combine these equations by premultiplying (9.14) by

4T

_{and solving (9 13) for u thus amving at a set of difference equations in}

standard form in x and A if we assume that Q and 4 are nonsingular These equa-tions are called Hamilton s equaequa-tions and their system matnx is the control Ham iltonian matrix:

[xl

[4t + 1'QI'MTQ1

-LAia+t I TQi lZ-T iLAi

The Haifliltonian for the control prOblem is, therefore,

e - _I

+ rQ;TT$:TQ1 -

_4,_TQ1 4,-T (9.29) Since theHarniltonian matrix is constant.. wecan solve (9.28) byfirst transforming

to a new state whichhas a diagonal system matrix, and from this solution easily

obtain the steady-state optimal control. In Appendix 9.1 we prove that the

eigen-values of this matrix are such that the reciprocal of every eigenvalue is also an eigenvalue. Therefore, half the roots of the characteristic equation must be inside

the unit circle and half must be outside. In this case, therefore, X can be diago-nalized to the form.4

E 0

=[

El'

where E is a diagonal matrix of the unstable roots _(IzI> 1) and E' is .a diagonal

matnx of the stable roots ((zI < I) W is obtained by the_{simlianty transforma}

tion

=

W'WW,

_(9.31)

where W is the matrix of eigenvectors OfW and can be written ifl block formas

w_lxt xc

IA,A0'

where

(9.30)

LA.

(9 27) is the matrix of eigenvectors associated with the eigenvalues (roots) outside the

unit circle and

(9.28)

_lxi

₌

_wixal

_Ix,

xoiFx

LA] LAJ LA, AOJLA

(9.32)

[xi

IA,

is the matrix of eigenvectors associated with the eigenvalues of _{which are} in-side the unit circle.

This same transformation matrix W can be used to transform x and A to the normal modes of the system, that is,

[:J

=

w'[jJ.

_(9.33)

where x and A are the normal modes. Conversely,_{we also have}

(9.34)

The solution to the coupled set of difference equations (9.28) can be simply stated in terms of the initial and final cofldtions and the normal modes, sincethe In rarecases, will haverepeatedroots and cannot be made diagonal by a change of variables. In

those cases a small ch.inge in Q, or Q will remove the problemor else we must compute the Jorddn

form for X. See Str.ing 0976). Reccnt method. Using the QZ algorithmare also possible.

258 MULTIVARIABLE AND OPTIMAL CONTROL

solution for the normal modeS is given by

11.1

fE_N olfxl

I.AJN 1. 0 _E"j1AJ0

To obtain the steady state, we let N go to infinity; therefore x, goes to zero

and in general A would go to infinity since each element of E is greater than one So we see that the only sensible solution for the steady-state (N x) case is for

A: = 0 and therefore A: = 0 for all k.5

From (9.34) and (935) with A: 0, we have

= X,x = X,Ex,

(9.36)

= A,x (9.37)

Therefore (9.36) leadS to

= EkXI¼,,. (9.38)

Thus, from (9.37) and (9.38),

A,XT'x = S.,x. (9.39)

Equation (9.39) is the same form as our assumption (9.16) for the sweep method so we conclude that S. given by (9 39) is the steady-state solution to (9.20) and that the control law for this system corresponding to with N mis

= Kx,

(9.40)

where, from (9.39) and (9.25),

(Q + rTsr)-'rs,.,'b.

(9.41)

Furthermore, from (927) the cost associated with using this control law is

= i4Sa.xo. (9.42)

The complete computational procedure is summarized5 below:

I Compute eigenvalues of the system matrix W defined by Eq (9 29)

Compute eigenvectors associatedwith the stable (Izi < I) eigenvalues and call

them

Ix,

LA,

Compute control gain IL from Eq. (9.41) with S. given by (9.39).

'From (9.34) we see that A is not zero then the state a will grow in time and the system will

beun-stable However if the system is controllable we know that a control exists which will make the system stable and give a finite value to Since we have the optimal control in (9 28) it must follow that the optimal system is stable and A' - 0 ii Q i is such that all the states affect5

9ti really necessary only to have a transformation that mskes X upper tnangular with all stable eigen values in the upper left half of the djagonal Sec Emani and Franklin (1979).

(9.35)

OPTIMAl. ESTIMATION 259

The stable eigenvalues from step I are the resulting system closed-loop roots

with constant gain K. from step 3 The reason for this is Before we let N go to infinity, the time-varying gain system was wntten in (9 29) as Imear constant

coefficient difference equations in x and A (9 29) The eigenvalues of this system

(of X) are equally valid for N small (time varying gain region)and N very large

(constant gain region) In the constant gain region it would be easy to uncouple

x from A by using (9.40), (9!41), and (9.10): and since eigenvalues remain un-changed for any linear combination of states, half of the eigenvalues of the W

system must be the eigenvahies of the closed-loop system

+ rue. (9.10) where

Uk = K,,Xa. (6.19)

We have already used the notion that the stable roots of W,. correspond tox. that

is (9 35) This is the only reasonable correspondence since the desired solution is that which minimizes the cost function (9 11) It therefore also follows that the

stable eigenvalues of the 4 system correspond to the Cigenvalues of the

closed loop system composed of(9 10) (9 40) and (9 42) We can also show that

the matrix X, of (9 34) and (9 36) is the matrix of eigenvectors of the optimal

teady-state closed-loop system,

An approach by Vaughn (1970) similar to the development in this section may be used to obtain a nonrecursive solution to the time varying case which could be considerably faster to solve than (9.21) through (9.25)

CON

MUL

stoorsignaal se..

CON

SEN

EQS.

SJN

CON

L

LOS..

MUL

DTNAMISCR POSITIONEE1NGSSYSTEM

noi: ruis . att: delen int:integratór

gai:vermeüigv. Con: constant, mul:vermenigv.

sum:soeer frq: sinusgenerätor.

sin:slnus COS.: cOsiflus...

NO

NOt

NOl

II I..

!ATT

Q

AT

ATT

INT

GAl

-10

INT

112

GAl

GAl

'I

y,

ir.

INT

tNT

INT

'3.

GAl

4$' schip en regeisysteem

X

lNT-

'C

GAi

I