Comments
Transcript
"Ingegneria dell`Automazione" - Volume 2
Manuale di Automazione by Cosimo Bettini, Timoti Lorenzetti, Alberto Guiggiani is licensed under a Creative Commons Attribuzione - Non commerciale - Condividi allo stesso modo 3.0 Unported License (http://creativecommons.org/licenses/by-ncsa/3.0/). di C.Bettini, A.Guiggiani e T.Lorenzetti 1 Indice I. Macchine ed azionamenti elettrici 1 1. Componenti 3 1.1. BJT . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.2. MOSFET . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 1.3. IGBT . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.4. TRIAC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 2. Convertitori DC-DC 2.1. Convertitore BUCK . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 9 2.1.1. Funzionamento continuo (CCM) . . . . . . . . . . . . . . . . . 10 2.1.2. Funzionamento discontinuo (DCM) . . . . . . . . . . . . . . . 13 2.1.3. Analisi del funzionamento limite . . . . . . . . . . . . . . . . . 17 2.1.4. Considerazioni . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 2.2. Convertitore Boost . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 2.2.1. Funzionamento continuo (CCM) . . . . . . . . . . . . . . . . . 20 2.2.2. Funzionamento limite . . . . . . . . . . . . . . . . . . . . . . . 24 2.2.3. Funzionamento discontinuo (DCM) . . . . . . . . . . . . . . . 25 2.2.4. Dimensionamento dei componenti . . . . . . . . . . . . . . . . 27 2.2.5. Considerazioni . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 2.3. Convertitore Cuk . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 2.3.1. Funzionamento continuo (CCM) . . . . . . . . . . . . . . . . . 30 2.3.2. Considerazioni . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 i Indice 2.4. Convertitore Sepic . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 2.4.1. Funzionamento continuo (CCM) . . . . . . . . . . . . . . . . . 34 2.4.2. Dimensionamento componenti . . . . . . . . . . . . . . . . . . 37 3. Convertitori isolati 41 3.1. Convertitore flyback . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 3.2. Convertitore forward . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 4. Analisi stazionaria dei convertitori 45 4.1. Esempio: Buck . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 5. Raddrizzatori 49 6. Azionamenti elettrici 65 6.1. Cenni sui trasformatori . . . . . . . . . . . . . . . . . . . . . . . . . . 65 6.1.1. Circuito reale . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 6.2. Motore asincrono . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 6.2.1. Motore asincrono monofase . . . . . . . . . . . . . . . . . . . 72 6.3. Motore sincrono . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72 II. Affidabilità e analisi dei rischi 77 7. Definizioni e concetti di base 79 7.1. Affidabilità . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 7.2. Disponibilità . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81 7.3. Prove . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87 7.3.1. Prove di affidabilità . . . . . . . . . . . . . . . . . . . . . . . . 87 8. Analisi dei rischi 89 8.1. Tecniche di analisi dei rischi . . . . . . . . . . . . . . . . . . . . . . . 91 9. Modelli e meccanismi di guasto 10.Safety e Security ii di C.Bettini, A.Guiggiani e T.Lorenzetti 95 101 Indice III. Sistemi Digitali 103 11.Conversione Analogico-Digitale (ADC) 105 11.1. Campionamento e quantizzazione . . . . . . . . . . . . . . . . . . . . 105 11.2. Convertitori . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 12.Conversione Digitale-Analogico (DAC) 119 12.1. Convertitori DAC e sintesi di segnali . . . . . . . . . . . . . . . . . . 123 13.Metastabilità e distribuzione di clock 127 13.1. Distribuzioni di clock . . . . . . . . . . . . . . . . . . . . . . . . . . . 131 14.Sistemi ad alta velocità 137 14.1. Terminazione di linee digitali . . . . . . . . . . . . . . . . . . . . . . . 141 15.Rumori e circuiti ad alta velocità 147 15.1. Switching-noise . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 147 15.2. Crosstalk . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 149 15.3. Problematiche di Layout . . . . . . . . . . . . . . . . . . . . . . . . . 152 IV. Relazione di macchine elettriche 15.4. Introduzione e Misure 157 . . . . . . . . . . . . . . . . . . . . . . . . . . 159 15.5. Obiettivi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 159 15.6. Dati di targa del motore trifase . . . . . . . . . . . . . . . . . . . . . 159 15.7. Misura della resistenza statorica . . . . . . . . . . . . . . . . . . . . . 159 15.8. Tabelle delle misure . . . . . . . . . . . . . . . . . . . . . . . . . . . . 160 15.9. Dati ricavati dalle misure . . . . . . . . . . . . . . . . . . . . . . . . . 160 16.Prova a vuoto 163 16.1. Tabelle delle misure . . . . . . . . . . . . . . . . . . . . . . . . . . . . 163 16.2. Dati ricavati dalla prova a vuoto . . . . . . . . . . . . . . . . . . . . 163 16.2.1. Separazione delle perdite nel ferro e meccaniche . . . . . . . . 165 16.3. Circuto equivalente a vuoto . . . . . . . . . . . . . . . . . . . . . . . 165 16.4. Grafici . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 166 di C.Bettini, A.Guiggiani e T.Lorenzetti iii Indice 16.5. Conclusioni . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 166 17.Prova in corto 169 17.1. Tabelle delle misure . . . . . . . . . . . . . . . . . . . . . . . . . . . . 169 17.2. Circuito equivalente in corto . . . . . . . . . . . . . . . . . . . . . . . 170 17.3. Parametri riportati alla temperatura convenzionale . . . . . . . . . . 171 17.4. Grafici . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172 17.5. Conclusioni . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172 18.Diagramma circolare 175 V. Laboratorio di automatica I 179 19.Controllo di un levitatore magnetico con PID e tecnica di taratura automatica IFT 181 20.Modello 183 20.1. Amplificatore di transconduttanza . . . . . . . . . . . . . . . . . . . . 183 20.2. Elettromagnete . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 184 20.3. Trasduttore di posizione . . . . . . . . . . . . . . . . . . . . . . . . . 186 20.4. Modello matematico del sistema . . . . . . . . . . . . . . . . . . . . . 187 20.5. Tempo di campionamento . . . . . . . . . . . . . . . . . . . . . . . . 188 21.Cenni sui regolatori PID 189 21.1. Il modello . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 189 21.2. Realizzazione dei regolatori PID . . . . . . . . . . . . . . . . . . . . 191 22.Iterative Feedback Tuning 193 22.1. Criteri di minimizzazione . . . . . . . . . . . . . . . . . . . . . . . . . 196 22.1.1. Calcolo dei gradienti . . . . . . . . . . . . . . . . . . . . . . . 197 22.2. Stima del gradiente 22.3. Convergenza . . . . . . . . . . . . . . . . . . . . . . . . . . . 200 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 201 22.4. Implementazione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 201 22.5. Calcolo dell’Hessiana di J . . . . . . . . . . . . . . . . . . . . . . . . 205 iv di C.Bettini, A.Guiggiani e T.Lorenzetti Indice 23.Costruzione dello schema Simulink 207 23.1. Schema Simulink . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 208 23.2. Primo Esperimento . . . . . . . . . . . . . . . . . . . . . . . . . . . . 208 23.3. Secondo Esperimento . . . . . . . . . . . . . . . . . . . . . . . . . . . 209 23.4. Algoritmo di aggiornamento . . . . . . . . . . . . . . . . . . . . . . . 210 23.5. Aggiornamento parametri . . . . . . . . . . . . . . . . . . . . . . . . 210 23.6. Segnale di controllo . . . . . . . . . . . . . . . . . . . . . . . . . . . . 211 24.Prove con riferimento a onda quadra 213 25.Prove con riferimento a dente di sega 219 26.Prove con riferimento sinusoidale 223 27.Robustezza del metodo 227 28.Conclusioni 231 VI. Laboratorio di automatica II 235 29.Controllo di un processo termico PT326 tramite PID con logica Fuzzy tarato attraverso algoritmi genetici 237 29.1. Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 237 29.2. Descrizione del Processo . . . . . . . . . . . . . . . . . . . . . . . . . 238 29.3. Identificazione del Processo . . . . . . . . . . . . . . . . . . . . . . . 240 30.PI+D con logica fuzzy 245 30.1. Discretizzazione del PID . . . . . . . . . . . . . . . . . . . . . . . . . 246 30.2. Controllore PI: . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 246 30.3. Controllore D: . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 247 30.4. Controllore PI+D: . . . . . . . . . . . . . . . . . . . . . . . . . . . . 249 31.Controllo Fuzzy 251 31.1. Fuzzificazione: . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 251 31.2. Regole d’inferenza: . . . . . . . . . . . . . . . . . . . . . . . . . . . . 253 di C.Bettini, A.Guiggiani e T.Lorenzetti v Indice 31.3. Defuzzificazione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 253 32.Algoritmi genetici 32.1. Introduzione agli algoritmi genetici . . . . . . . 32.2. Definizione del problema . . . . . . . . . . . . . 32.3. Operazioni genetiche . . . . . . . . . . . . . . . 32.4. Applicazioni per l’ottimizzazione del controllore 32.5. Implementazione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33.Prove 259 . 259 . 260 . 260 . 262 . 263 267 VII.Elaborati di stima e identificazione 283 34.Applicazioni del filtro di Benes 34.1. Introduzione . . . . . . . . . 34.2. Modello del sistema . . . . . 34.3. Filtro di Benes . . . . . . . 34.4. Filtraggio EKF . . . . . . . 34.5. Filtraggio PF . . . . . . . . 34.6. Prove Monte Carlo . . . . . . . . . . . . . . . . . 285 . 285 . 288 . 289 . 297 . 302 . 318 35.Stima dello stato in presenza di 35.1. Modellazione . . . . . . . . 35.2. Algoritmo EKF . . . . . . . 35.3. Algoritmo UKF . . . . . . . 35.4. Conclusioni . . . . . . . . . osservazioni binarie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 323 . 323 . 326 . 331 . 337 vi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . di C.Bettini, A.Guiggiani e T.Lorenzetti . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Parte I. Macchine ed azionamenti elettrici 1 1. Componenti Iniziamo con l’analisi dei vari componenti che vengono utilizzati nei circuiti elettronici di potenza per la costruzione dei circuiti per l’azionamento delle macchine elettriche. Si noti che in componenti elettronici attivi quali i transistor nell’elettronica di potenza vengono utilizzati nel loro funzionamento non lineare, ovvero come interruttori. 1.1. BJT Il BJT (bipolar junction transistor) è un particolare transistor in cui tre strati di materiale semiconduttore drogato in cui lo strato centrale detto base è drogato in maniera opposta agli altri due (detti collettore ed emettitore) in modo da formare due giunzioni p-n. Si ottengono così due configurazioni di BJT: BJT p-n-p, BJT n-p-n. Nell’analisi seguente si considera in BJT n-p-n (Figura 1.1). Il BJT lavora come interruttore se è in condizione di interdizione/saturazione, il comando per poter pilotare il BJT nelle due fasi è la corrente di base IB . Quando il transistor è spento (interdetto IB = 0) la tensione VCE è pari a VCC applicata e la corrente IC = 0 ne consegue che la potenza dissipata è pari a POF F = VCE IC = 0. Quando il transistor è acceso (saturazione IB > 0) la tensione VCE è circa zero (VCE ' 0.3) e scorre una corrente IC in base al carico pilotato ne consegue che la potenza dissipata è pari a PON = VCE IC 0, introducendo problemi di surriscaldamento del componente. Inoltre la frequenza di commutazione del BJT è relativamente bassa (ORDINE). 3 1.2 MOSFET Figura 1.1.: BJT Il comportamento nella fase di ON del BJT è modellizzabile come un condensatore, per cui per passare alla fase di OFF si deve “aspettare” che la capacità si scarichi. Ecco perchè il BJT risente di tempi di commutazione relativamente alti (ha un tempo di spegnimento alto). 1.2. MOSFET Il MOSFET (metal oxide semiconductor field effect transistor) è un particolare transistor composto da un substrato di materiale semiconduttore drogato a cui sono applicati tre terminali: source, gate, drain (Figura 1.2). In base all’applicazione di una tensione al gate permette di controllare il passaggio di corrente tra il source e il drain. Quando il transistor è spento (interdetto VGS ≤ VT 0 con VT 0 tensione di soglia) non scorre corrente tra surce e drain (ISD = 0) e quindi anche in questo caso la potenza in condizione di interdizione è pari a zero. Quando il transistor è acceso (saturazione VDS = VGS − VT 0 ) la corrente IDS è libera di scorrere. Il comportamento nella fase di ON è modellizzabile come una resistenza 2 (RDSON ) ne consegue che la potenza dissipata è PON = RDSON IDS . 4 di C.Bettini, A.Guiggiani e T.Lorenzetti 1.3 IGBT Figura 1.2.: MOSFET Il MOSFET risulta essere sicuramente più veloce rispetto al BJT nei tempi di commutazione, infatti in generale generale i BJT vengono utilizzati nei circuiti di potenza (dissipano meno), mentre i MOSFET vengono utilizzati nei circuiti di pilotaggio (dinamica veloce). 1.3. IGBT L’IGBT (insulate gate bipolar transistor) è un transistor che fonde le caratteristiche di un BJT e di un MOSFET: nell stadio di ingresso ha le caratteristiche del MOSFET e nello stadio di uscita quelle del BJT (Figura 1.3). Figura 1.3.: IGBT di C.Bettini, A.Guiggiani e T.Lorenzetti 5 1.4 TRIAC 1.4. TRIAC Il TRIAC è un componete semiconduttore il cui funzionamento può essere considerato come un diodo controllato, infatti è un dispositivo con tre terminali (Figura 1.4): due anodi tra i quali è possibile fare scorrere la corrente ed il gate che funge da ingresso di controllo (attivazione). Infatti ciascun “diodo” conduce nel semiperiodo dell’onda elettrica i cui è polarizzato direttamente in accordo se viene applicata una corrente superiore alla soglia di sensibilità sul gate. Figura 1.4.: TRIAC 6 di C.Bettini, A.Guiggiani e T.Lorenzetti 2. Convertitori DC-DC I convertitori DC-DC sono dispositivi molto utilizzati nell’elettronica di potenza e servono a convertire un livello di tensione in continua in un altro livello. I convertitori si dividono in due grandi famiglie: convertitori isolati e convertitori non isolati. I convertirori isolati mantengono isolati gli stadi di ingresso e di uscita del convertitore grazie ad un trasformatore di isolamento galvanico raggiungendo livelli di sicurezza migliori in caso di guasto. Inoltre riescono ad ottenere rapporti di conversione più alti rispetto ai convertirori non isolati. Tutti i convertitori DC-DC si basano sul circuito chopper caretterizzato dal principio di funzionamento della commutazione. Si consideri il circuito in Figura 2.1 formato da un interruttore ed un diodo con applicata in ingresso una tensione continua. Figura 2.1.: Chopper In uscita dal convertitore l’andamento della tensione risulta come in Figura 2.2, il valor medio della tesione è dato da: 7 Convertitori DC-DC VOU T = 1 T T ˆON vout (t)dt = T1 VI TON 0 Si definisce duty cicle il rapporto D = TON , T con T periodo di commutazione, quindi: VOU T = DVI Figura 2.2.: Andamento della tensione di uscita dal Chopper Il valor medio della tensione di uscita risulata essere proporzionale al duty cicle, ovvero è possibile modificare il valore della tensione di uscita modificando i tempi di apertura o di chiusura dell’interruttore. Questo tipo di conversione, basata sullo switching (commutazione), risulta essere molto più efficiente rispetto ad un semplice partitore che è un convertiore statico, dissipativo e in cui non è possibile inalzare il livello di tensione. In generale i convertitori che utilizzano questo principio di funzionamento studiando varie metodologie (filtri) per poter estrapolare la componente in continua (valor medio) della tensione in uscita dal circuito chopper. Per poter analizzare il funzionamento dei convertitori si fa riferimento alle seguenti tre ipotesi: 1. funzionamento stazionario 2. alimentazione continua ideale 8 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.1 Convertitore BUCK 3. carico resistivo o homico-induttivo Inoltre durante l’analisi si fa riferimento anche al modo di funzionamento, ovvero si considera il comportamento della corrente sull’induttore (o la tensione sul condensatore) del filtro: CCM si dice che il convertitore è in funzionamento continuo (continuos conduction mode) se la corrente sull’induttore (o la tensione sul condensatore) non cambia di polarità durante il ciclo di commutazione. In questo caso la tensione di uscita non dipende dai componenti circuitali. DCM si dice che il convertitore è in funzionamento discontinuo (discontinuos conduction mode) se la corrente sull’induttore (o la tensione sul condensatore) risulta essere zero per un certo tempo durante il ciclo di commutazione. In questo caso la tensione di uscita dipende dai componenti circuitali. 2.1. Convertitore BUCK Il convertitore Buck appartiene alla classe di convertitori non isolati. Esso permette di ottenere un livello di tensione più basso di tensione in uscita rispetto a quello in ingresso. Infatti viene detto anche convertitore step-down, ovvero abbassatore. Lo schema di principio lo si può notare in Figura 2.3. Figura 2.3.: Buck di C.Bettini, A.Guiggiani e T.Lorenzetti 9 2.1 Convertitore BUCK 2.1.1. Funzionamento continuo (CCM) Fatte le precedenti considerazioni possiamo ricavare i circuiti equivalenti nelle fasi di ON e di OFF rappresentati in figura Figura 2.4, mentre per quanto riguarda l’andamento della corrente e la tensione sull’induttore è ragionevole ipotizzare un andamento del tipo illustrato in figura Figura 2.5. Figura 2.4.: Circuiti equivalenti del Buck in condizione rispettivamente di ON e di OFF 10 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.1 Convertitore BUCK Figura 2.5.: Curve caratteristiche di un convertitore Buck in CCM 2.1.1.1. Fase di ON: interruttore chiuso Come si può notare in figura, durante questa fase la tensione sull’induttore è pari a: VL = VIN − Vout mentre la corrente IL cresce linearmente. Essendo l’interruttore chiuso, il diodo viene polarizzato inversamente e non vi è quindi circolazione di corrente su di esso. Il valore della corrente (IL ) è dato dalla relazione: VL = L · diL dt di C.Bettini, A.Guiggiani e T.Lorenzetti 11 2.1 Convertitore BUCK Ne deriva così, invertendo la relazione ed integrando, l’incremento (∆IL ) di corrente durante la fase di ON è pari a: iˆ max ∆ILon = ˆton diL = 0 imin VL (VIN − Vout ) ton dt = L L 2.1.1.2. Fase di OFF: interruttore aperto Nella fase di OFF, invece, la variazione di corrente è pari a: iˆ max ∆ILof f = ˆtof f diL = 0 imin VL Vout · tof f dt = − L L Assumendo che il convertitore lavori in regime stazionario, ne consegue che l’energia immagazzinata in ciascun componente alla fine del ciclo di commutazione debba essere uguale a quella di inizio ciclo. Questo significa che il valore della corrente iL è lo stesso per t = 0 e t = T (dove T è il periodo di commutazione). Possiamo quindi scrivere che: ∆ILon − ∆ILof f = (VIN − Vout ) ton Vout · tof f − =0 L L da cui, dato: D= ton , T = ton + tof f T si ottiene: 12 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.1 Convertitore BUCK Vout = D · VIN Da questa equazione si può osservare come la tensione di uscita del convertitore vari linearmente con il duty cycle (il quale non può risultare maggiore dell’unità). Per questo motivo questo convertitore è noto anche col nome di convertitore step-down. 2.1.2. Funzionamento discontinuo (DCM) Lavorando in DCM si assiste, di fatto, alla creazione un terzo stato in cui l’interruttore è aperto e, al contempo, il diodo è interdetto (quindi la corrente IL = 0). Figura 2.6.: Circuito equivalente del Buck in condizione del terzo stato L’andamento della corrente e della tensione sull’induttore sarà come quella in figura: di C.Bettini, A.Guiggiani e T.Lorenzetti 13 2.1 Convertitore BUCK Figura 2.7.: Curve caratteristiche di un convertitore Buck in DCM in cui è possibile evidenziare i tre stati di funzionamento: • [0 , DTS ] Fase di ON; • [DTS , D1 TS ] Fase di OFF; • [D1 TS , TS ] Fase di OFF dell’interruttore e diodo. Possiamo notare che, a parità di energia trasferita, la corrente IM AX nel caso di funzionamento DCM è più alta rispetto al caso CCM (dovendo essere le aree descritte dall’andamento delle correnti IL equivalenti). Nel caso discontinuo, le espressioni delle variazioni di corrente sono ricavabili dalle curve caratteristiche riportate in figura Figura 2.7. Più precisamente: ∆IL+ 14 = 1 D (Vi − Vo ) TS L di C.Bettini, A.Guiggiani e T.Lorenzetti 2.1 Convertitore BUCK ∆IL− = 1 Vo (D1 − D) TS L Il valor medio delle corrente risulta così: 1 1 1 = IM AX D1 IL = IM AX D1 TS · 2 TS 2 La corrente sull’induttore è nulla all’inizio e cresce durante il periodo ton . Come si può notare in figura Figura 2.7, vale: IM AX = Vin − Vo DT L Sostituendo e assumendo che la corrente media sull’induttore è pari alla corrente media in uscita: Io = IL = (Vin − Vo )DT D1 2L Come esposto in precedenza, deve valere ∆IL+ = ∆IL− (figura Figura 2.7) e si ottiene: D Vo (Vin − Vo ) T = (D1 − D) T L L da cui: DVin = Vo D1 di C.Bettini, A.Guiggiani e T.Lorenzetti 15 2.1 Convertitore BUCK Definendo mV come il guadagno ingresso-uscita, si ottiene: mv = Vo D = Vin D1 Sappiamo inoltre che nel buck I0 = IL ; quindi: D1 = 2I0 ∆I+ da cui è possibile esplicitare D1 come: D1 = 2I0 S (VI − V0 ) D·T L Quindi, riprendendo l’equazione del guadagno mV e sostituendo l’espressione di D1 appena ricavata, si ottiene: mV = D D (VI − V0 ) DTS = = D1 2I0 L D2 1 mV −1 R 2Lf da cui, esplicitando D: s D = mv 2Lf 1 R 1 − mV Si osservi infine che la relazione tra tensione di ingresso e di uscita (mv ) non è più dipendente solo dal duty cycle (come nel caso CCM) ma anche dai componenti del circuito, rendendo così la caratterizzazione più complessa ma, al contempo, più varia. 16 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.1 Convertitore BUCK 2.1.3. Analisi del funzionamento limite Per concludere l’analisi del convertitore Buck andiamo adesso ad analizzare cosa accade al confine tra funzionamento CCM e DCM. E’ possibile identificare tale zona come quella in cui vale la relazione: D1 T = T ⇒ D1 = 1 Al contempo: Iolimite = D1 IM AX Vin − Vo IM AX = = DT 2 2 2L che, considerando la relazione: Vo = DVin implica: Iolimite = Vin (1 − D) DT 2L Introduciamo adesso la tensione e la corrente normalizzate definendole come: |Vo | = Vo L , |Io | = Io Vin Vin T Ne consegue che in CCM: di C.Bettini, A.Guiggiani e T.Lorenzetti 17 2.1 Convertitore BUCK |Vo | = D mentre in DCM: |Vo | = 1 2LIo D2 Vin T +1 = D2 = 2 |Io | + D2 +1 1 2|Io | D2 con la corrente limite Iolimite pari a: Iolimite = Vin Io D(1 − D)T = D(1 − D) 2L 2 |Io | Quindi, nel punto limite tra funzionamento DCM e CCM si verifica: Iolimite = Io da cui: (1 − D)D =1 2 |Io | per comprendere meglio tali relazioni si fa riferimento alla figura Figura 2.8. 18 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.1 Convertitore BUCK Figura 2.8.: Relazione tra corrente e tensione normalizzate in DCM e CCM. 2.1.4. Considerazioni Se i componenti fossero ideali si avrebbe: PIN = POU T → VIN IL = VOU T IOU T → VOU T IL = IOU = D. In verità la relazione statica VOU T = DVIN è valida per D = VIN T (0.3 ÷ 0.7) per effetto delle perdite dovute ai componenti non ideali. di C.Bettini, A.Guiggiani e T.Lorenzetti 19 2.2 Convertitore Boost Figura 2.9.: Funzione di trasferimento del Buck 2.2. Convertitore Boost Il convertitore Boost appartiene alla classe di convertitori non isolati. Esso permette di ottenere un livello di tensione più alto uscita rispetto a quello in ingresso. Infatti viene detto anche convertitore step-up, ovvero alzatore. Lo schema di principio lo si può notare in Figura 2.10. Figura 2.10.: Boost 2.2.1. Funzionamento continuo (CCM) Andiamo ad analizzare i circuiti equivalenti nelle fasi di ON e di OFF rappresentati in Figura 2.11, mentre per quanto riguarda l’andamento della corrente e della tensione 20 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.2 Convertitore Boost sull’induttore si ipotizza come in Figura 2.12. Figura 2.11.: Circuiti equivalenti nelle fasi ON/OFF del Boost di C.Bettini, A.Guiggiani e T.Lorenzetti 21 2.2 Convertitore Boost Figura 2.12.: Tensione e corrente sull’induttore 2.2.1.1. Fase di ON: interruttore chiuso Durante questa fase la tensione sull’induttore è pari a: vL = VIN quindi la corrente sull’induttore (iL ) è data dalla relazione: diL vL = L dt cresce linearmente con pendenza VIN . L La corrente sul condensatore invece è data dalla relazione: 22 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.2 Convertitore Boost iC = − VOU T R 2.2.1.2. Fase di OFF: interruttore aperto Le relazioni che si ottengono nella fase di OFF sono: vL = VIN − VOU T = L iC = iL − diL dt VOU T R Come si può notare dalla Figura 2.12 per essere nel funzionamento CCM deve valere iL (0) = iL (T ) > 0, quindi considerando le equazioni: iL (DT ) − iL (0) = VIN DT L iL (DT ) − iL (0) = −( VIN − VOU T )(1 − D)T L si ottiene la relazione tra la tensione di ingresso e di uscita caratteristica del convertitore Boost: VOU T = VIN 1−D di C.Bettini, A.Guiggiani e T.Lorenzetti 23 2.2 Convertitore Boost 2.2.2. Funzionamento limite Per prima cosa si consideri il fatto che se i componenti fossero ideali la potenza di ingresso e di uscita sarebbero uguali: PIN = POU T 2 IM IN + IM AX VOU T VIN = 2 R inoltre si nota che: ∆iL = IM AX − IM IN = VIN DT L quindi è possibile ricavarsi le seguenti due espressioni: IM AX = VIN R(1−D)2 + VIN DT 2L IM IN = VIN R(1−D)2 − VIN DT 2L (2.1) La condizione di funzionamento limite si ha quando IM IN = 0. Quindi la corrente limite in ingresso al conertitore è data da: IIN,limite = ∆iL VIN DT = 2 2L mentre la corrente limite in uscita è data da: 2 IOU T,limite = (1 − D) IIN,limite 24 VIN DT (1 − D)2 = 2L di C.Bettini, A.Guiggiani e T.Lorenzetti 2.2 Convertitore Boost dato che la corrente “passa” solo nella fase di OFF dell’interruttore. 2.2.3. Funzionamento discontinuo (DCM) Andiamo ad analizzare i circuiti equivalenti nelle fasi di ON e di OFF rappresentati in Figura 2.11. Nel caso discontinuo (DCM) l’andamento di tensione e corrente corrisponde a quello in Figura 2.13: la corrente si annulla in D1 T e resta zero fino a T , quando riprende a salire. Dato che in CCM diodo e interruttore commutano per Imin 6= 0, ciò porta alla presenza di spike di corrente dovuti alla commutazione; in questo caso, invece,accensione e spegnimento di interruttore e diodo avvengono in I = 0. Di contro, il picco di corrente IM AX è maggiore nel caso DCM rispetto al CCM (e ciò va tenuto di conto nella fase di progetto). Andiamo adesso ad analizzare la relazione tra D, corrente e tensione del boos nel caso discontinuo. Si ottiene: VIN D = (VOU T − VIN )(D1 − D) da cui si ottiene: D1 VOU T = VIN D1 − D da cui deriva: 1 1 VIN D IL = D1 IM AX = D1 2 2 Lf e, di conseguenza: D1 = 2 Lf IL Lf D1 =2 · IOU T VIN D VIN D1 − D con f la frequenza di funzionamento del convertitore. Si ottiene inoltre: 1 D1 − D IOU T = ID = IM AX (D1 − D) = IL 2 D di C.Bettini, A.Guiggiani e T.Lorenzetti 25 2.2 Convertitore Boost e di conseguenza si definisce: ·IOU T D + 2 L·f VOU T D·VIN = Mv = ·IOU T VIN 2 L·f D·VIN −1 che dipende interamente da D. Considerando D · VIN = VOU T ed R = VOU T IOU T si ottiene la formula finale: D2 = 2 L·f Mv (Mv − 1) R che esprime l’andamento di D in funzione di Mv . Introducendo la notazione: −1 • V̄ è la tensione normalizzata, definita dal rapporto VOU T VIN e corrispondente al guadagno in tensione del convertitore; • I¯ è la corrente normalizzata, definita dalla relazione I¯ = T VLIN IOU T dove T VLIN è pari all’incremento massimo della corrente dell’induttore durante un ciclo, cioè l’incremento della corrente dell’induttore con un duty cycle D = 1, quindi in regime stazionario del convertitore (ciò significa che è uguale a 0 per corrente di uscita nulla e uguale ad 1 per la massima corrente che il convertitore può fornire). Usando queste notazioni avremo che: • in modo continuo: V̄ = 1 1−D • in modo discontinuo: V̄ = 1 + VIN D2 T D2 =1+ ¯ 2LIOU T 2I • nel punto limite tra modo continuo e discontinuo si ha che la corrente limite tra modo continuo e discontinuo è data da: Ilim = 26 VIN T Ilim D(1 − D) = ¯ D(1 − D) 2L 2I di C.Bettini, A.Guiggiani e T.Lorenzetti 2.2 Convertitore Boost perciò il punto limite tra modo continuo e discontinuo è dato da: 1 D(1 − D) = 1 2I¯ Queste espressioni sono state tracciate in Figura 2.14. Le differenza di comportamento tra modo continuo e discontinuo può essere chiaramente osservata. Ciò è molto importante dal punto di vista di un circuito di controllo. Figura 2.13.: Funzionamento DCM del boost 2.2.4. Dimensionamento dei componenti Considerando l’equazione Equazione 2.1 e che il funzionamento limite si ha quando IM IN = 0 ci si può ricavare l’espressione dell’induttanza critica: Lcrit = RT D(1 − D)2 2 di C.Bettini, A.Guiggiani e T.Lorenzetti 27 2.2 Convertitore Boost Figura 2.14.: Funzionamento DCM del boost quindi se si vuole che il convertitore “funzioni” in CCM o DCM si dovrà scegliere un’induttanza maggiore o minore a quella critica. Per quanto riguarda il condensatore si consideri il ripple di tensione su di esso in accordo a gli schemi in Figura 2.11 e considerando il funzionamento CCM (vedi Figura 2.15). Figura 2.15.: Ripple di tensione sul condensatore 28 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.2 Convertitore Boost Si ha: ∆Vripple = − − VOU T RC DT quindi: C= VOU T DT R∆Vripple In questo caso quinidi ∆Vripple è un parametro di progetto. 2.2.5. Considerazioni Se in componenti fossero ideali per D → 1 la tensione VOU T → ∞ grazie alla VIN . In verità questa relazione in generale vale fino a relazione statica VOU T = 1−D DM AX = 0, 9 a causa della non idealità dei componenti (vedi Figura 2.16). Figura 2.16.: Comportamento reale del Boost di C.Bettini, A.Guiggiani e T.Lorenzetti 29 2.3 Convertitore Cuk 2.3. Convertitore Cuk Il convertitore Cuk appartiene alla classe di convertitori non isolati. Esso permette di ottenere un livello di tensione più alto o più alto in uscita rispetto a quello in ingresso. Infatti lo schema di principio come si può notare in FIGURACUK è una combinazione tra gli schemi del Boost e del Buck. Figura 2.17.: Buck 2.3.1. Funzionamento continuo (CCM) Andiamo ad analizzare i circuiti equivalenti nelle fasi di ON e di OFF rappresentati in Figura 2.18. 30 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.3 Convertitore Cuk Figura 2.18.: Circuiti equivalenti nelle fasi ON/OFF del Cuk 2.3.1.1. Fase di ON: interruttore chiuso Durante questa fase le tensioni sugli induttori L1 e L2 è pari a: vL1 = VIN vL2 = −VOU T − VC1 2.3.1.2. Fase di OFF: interruttore aperto Le relazioni che si ottengono nella fase di OFF sono: vL1 = VIN − VC1 di C.Bettini, A.Guiggiani e T.Lorenzetti 31 2.3 Convertitore Cuk vL2 = −VOU T Figura 2.19.: Andamento della tensione sugli induttori Come si può notare dalla Figura 2.19 per essere nel funzionamento CCM le due aree si devono equivalere, quindi: VIN D = (VIN − VC1 ) (1 − D)2 (VC1 − VOU T ) D = VOU T (1 − D) 32 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.4 Convertitore Sepic si ottiene così la relazione tra la tensione di ingresso e di uscita caratteristica del convertitore Cuk: VOU T = − D VIN 1−D come si può notare il convertitore Cuk risulta essere invertente. 2.3.2. Considerazioni Il convertitore Cuk ha il vantaggio di avere ripple di corrente di ingreresso e di usicta contenuti e presenta un naturale isolamento capacitivo contro i guasti dell’interruttore. Di contro il condensatore C1 ha un’elevata corrente di ripple e quinidi è presente un’elevata tensione e corrente sull’interuttore . 2.4. Convertitore Sepic Il convertitore Sepic appartiene alla classe di convertitori non isolati. Esso permette di ottenere un livello di tensione più alto o più alto in uscita rispetto a quello in ingresso. Rispetto al convertitore Cuk non è invertente. Figura 2.20.: Buck di C.Bettini, A.Guiggiani e T.Lorenzetti 33 2.4 Convertitore Sepic 2.4.1. Funzionamento continuo (CCM) Andiamo ad analizzare i circuiti equivalenti nelle fasi di ON e di OFF rappresentati in Figura 2.21, mentre per quanto riguarda l’andamento della tensione sugli induttori si ipotizza che sia come nella Figura 2.22 Figura 2.21.: Circuiti equivalenti nelle fasi ON/OFF del Sepic 2.4.1.1. Fase di ON: interruttore chiuso Durante questa fase le tensioni sugli induttori L1 e L2 sono pari a: vL1 = VIN vL2 = VC1 34 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.4 Convertitore Sepic Durante questa fase le correnti sui condensatori C1 e C2 sono pari a: iC1 = −iL2 iC2 = − VOU T R 2.4.1.2. Fase di OFF: interruttore aperto Le relazioni che si ottengono nella fase di OFF sono: vL1 = VIN − VC1 − VOU T vL2 = −VOU T iC1 = −iL1 iC2 = − VOU T + iL1 + iL2 R di C.Bettini, A.Guiggiani e T.Lorenzetti 35 2.4 Convertitore Sepic Figura 2.22.: Andamento della tensione sugli induttori Come si può notare dalle Figura 2.22 per essere nel funzionamento CCM le due aree si devono equivalere, quindi: VIN DT VC DT 1 = (VIN − VC1 − VOU T ) (1 − D) = −VOU T (1 − D)T −i DT = iL1 (1 − D)T L2 − VOU T DT = i + i − VOU T (1 − D)T R L1 L2 R risolvendo il sistema si ottiene la relazione tra la tensione di ingresso e di uscita caratteristica del convertitore Sepik: 36 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.4 Convertitore Sepic VOU T = D VIN 1−D come si può notare il convertitore Cuk risulta essere invertente. 2.4.2. Dimensionamento componenti Per considerare il dimensionamento delle induttanze ci si va sempre a ricavare l’induttanza critica, ovvero ci si ricava l’induttanza nel caso di funzionamento limite quindi: ∆iL = IL1 M AX − IL1 M IN = IL M AX 1 = IL1 + VIN DT 2L1 IL1 M IN = IL1 + VIN DT 2L1 VIN DT L1 ponendo IL1 M IN = 0 si ricava: L1crit = R(1 − D)2 T 2D allo stesso modo per L2 : ∆iL = IL2 M AX − IL2 M IN = VC1 DT L1 di C.Bettini, A.Guiggiani e T.Lorenzetti 37 2.4 Convertitore Sepic IL M AX 2 = IL2 + VIN DT 2L2 IL2 M IN = IL2 + VIN DT 2L2 ponendo IL2 M IN = 0 si ricava: L2crit = R(1 − D) T 2 Per quanto riguarda i condensatori si consideri il ripple di tensione su di essi, in accordo a gli schemi in Figura 2.21 e considerando il funzionamento CCM (si osservi Figura 2.23). Figura 2.23.: Ripple di tensione sui condensatori 38 di C.Bettini, A.Guiggiani e T.Lorenzetti 2.4 Convertitore Sepic ∆VC1 ripple = IL2 DT C1 ∆VC2 ripple = VOU T DT RC2 quindi: C1 = C2 = IL2 ∆VC1 ripple DT VOU T DT R∆VC2 ripple di C.Bettini, A.Guiggiani e T.Lorenzetti 39 3. Convertitori isolati Analizziamo adesso i convertitori isolati. Questi sono ottenuti ponendo un trasformatore all’interno del convertitore in modo tale da ottenere una separazione galvanica delle cariche. I convertitori isolati si dividono in due categorie molto differenti: convertitori simmetrici e convertitori asimmetrici. La simmetria è riferita alla posizione del punto di lavoro magnetico del trasformatore (il piano B − H) rispetto all’origine degli assi; sono quindi simmetrici quelli a due quadranti e asimmetrici quelli ad un quadrante. 3.1. Convertitore flyback Il convertitore flyback è un derivato dal buck-boost ed è caratterizzato dalla presenza di un singolo interruttore (mentre, ad esempio, nei convertitori half-bridge e push-pull sono presenti due interruttori); questo convertitore è molto usato nell’alimentazione dei personal computer in quanto adatto ad applicazioni a basso costo e bassa potenza. Se applichiamo il concetto di isolamento galvanico ad un convertitore buck-boost si ottiene il risultato in Figura 3.1: il problema è che con questa configurazione è presente una induttanza di dispersione Ld che serve a tener conto del fatto che non tutto passa da primario a secondario; quando l’interruttore apre, l’induttanza Ld non può (per come è stata definita) trasferire energia al secondario, provocando delle sovratensioni che vanno ad agire sulll’interruttore (che può addirittura bruciarsi). Per sopperire a questa problematica è possibile introdurre dei diodi di ricircolo od un apposito circuito di snubber (in Figura 3.2), considerando però che con queste confi- 41 3.2 Convertitore forward gurazioni l’energia in eccesso viene dissipata su una resistenza e quindi il rendimento del convertitore si abbassa (e non di poco). Si pensa adesso di scaricare l’energia in eccesso sul primario così da rigenerarla e mantenere comunque l’isolamento. Per fare ciò si introduce un avvolgimento di recupero (in Figura 3.3), naturalmente con sezione minore di quella del primario e tale che, se ben modellato, l’energia in eccesso venga re-immessa sul primario. In questo modo sono sì presenti delle sovratensioni, però è possibile essere consapevoli della loro entità e limitarle modellando l’avvolgimento di recupero. Figura 3.1.: Convertitore flyback isolato Figura 3.2.: Snubber per converitore flyback isolato 3.2. Convertitore forward Il convertitore forward è un derivato dal buck e può essere immaginato come in Figura 3.4. Il problema di questa configurazione è che, a causa dell’induttanza 42 di C.Bettini, A.Guiggiani e T.Lorenzetti 3.2 Convertitore forward Figura 3.3.: Convertitore flyback isolato con avvolgimento di recupero di magnetizzazione Lm , la corrente circolante su di essa va costantemente a salire (quando l’interruttore è chiuso, mentre è costante a interruttore aperto in quanto VLm = 0) fino a saturare raggiungendo valori molto elevati (valori che prima della fine fanno saltare l’interruttore). È necessario utilizzare un circuito tale che la tensione sull’induttanza Lm sia negativa quando l’interruttore è chiuso (così da tenere limitata la tensione iLm ); per fare ciò si può utilizzare la logica in Figura 3.5 andando ad accoppiare l’induttanza del secondario con un terzo ramo del trasformatore così da avere una tensione a interruttore aperto del tipo: VLm = − N3 Vout N Dato che l’andamento della tensione sull’induttanza è come quello in Figura 3.6, per l’uguaglianza delle aree si ottiene: N Vout D = Vin (1 − D) N3 che ci porta alla condizione: N3 Vin D > N Vout (1 − D) che deve essere soddisfatta affinchè quanto detto abbia senso. Si noti inoltre il diodo sul ramo aggiunto necessario a fare sì che l’espediente introdotto funzioni solo quando l’interruttore è aperto (quando è chiuso perderebbe di utilità). di C.Bettini, A.Guiggiani e T.Lorenzetti 43 3.2 Convertitore forward Figura 3.4.: Convertitore forward isolato Figura 3.5.: Convertitore forward isolato Figura 3.6.: Andamento della tensione su Lm in un convertitore forward isolato 44 di C.Bettini, A.Guiggiani e T.Lorenzetti 4. Analisi stazionaria dei convertitori Nell’analisi stazionaria dei convertitori è difficile capire come trattare l’interruttore che risulta essere un dispositivo tempo variante. Vediamo come analizzare tale componente considerato in Figura 4.1. Figura 4.1.: Interruttore La relazione statica dell’interruttore è data da VLD = DVSD . Per analizzare il comportamento tempo variante dell’interruttore, si fa l’ipotesi di considerare la variazione temporale come una perturbazione (vedi Figura 4.2). Quindi la relazione statica diventa: (VLD + vld ) = (D + d) (VSD + vsd ) sviluppando: 45 Analisi stazionaria dei convertitori Figura 4.2.: Approssimazione della variazione del duty cicle come perturbazione (VLD + vld ) = DVSD + Dvsd + dVSD + dvsd trascurando i termini di “secondo ordine” (dvsd = 0) possiamo scomporre la tensione in uscita all’interruttore nella parte statica (DVSD ) e nella parte alle variazioni (Dvsd + dVSD ). Lo stesso ragionamento viene applicato alla corrente considerando la relazione statica IS = DIL . (IS + is ) = (D + d) (IL + il ) = DIL + Dil + dIL + dil trascurando il termine dil si ottiene il circuito equivalente dell’interruttore (vedi Figura 4.3) 46 di C.Bettini, A.Guiggiani e T.Lorenzetti Analisi stazionaria dei convertitori Figura 4.3.: Modello della cella a commutazione Si può mettere in evidenza il comportamento stazionario (Figura 4.4) e alle variazioni (Figura 4.5) dell’interruttore. Figura 4.4.: Modello stazionario Figura 4.5.: Modello alle variazioni di C.Bettini, A.Guiggiani e T.Lorenzetti 47 4.1 Esempio: Buck 4.1. Esempio: Buck Si introduce ora l’analisi stazionaria e alle variazioni del covertitore Buck. Si ricorda che in regime stazionario l’induttore si comporta come un circuito chiuso mentre il condensatore come un circuito aperto, quindi in accordo a quanto detto precedentemente possiamo ricavarci il modello equivalente del Buck in regime stazionario (Figura 4.6). Figura 4.6.: Analisi stazionaria del Buck VOU T = DVSD = DVIN Alle variazioni invece si ottiene il modello rappresentato in Figura 4.7. Figura 4.7.: Analisi stazionaria del Buck 48 di C.Bettini, A.Guiggiani e T.Lorenzetti 5. Raddrizzatori In elettrotecnica ed elettronica un raddrizzatore (rectifier in inglese) è un dispositivo che serve a raddrizzare un segnale bipolare trasformandolo in un segnale unipolare . Il raddrizzatore, collegato ad altri componenti, è usato quindi per trasformare la corrente alternata (AC) in corrente continua (DC). Nel fare ciò si cerca di massimizzare il rendimento e, al contempo, ottenere una uscita il più stabile possibile. Si noti che nelle norme la tensione DC è tale per cui il suo ripple non supera il 10% Raddizzatore a singola semionda il segnale d’ingresso, sinusoidale, viene applicato ad un diodo in serie alla resistenza di carico (Figura 5.4.a e Figura 5.5.a). Se il catodo è rivolto verso il carico, esso consente il passaggio delle sole semionde positive lasciando uno spazio tra semionda positiva e semionda positiva perché il raddrizzatore a semionda non trasforma le onde negative in positive. Notiamo che in un raddrizzatore del genere, data la tensione di ingresso: vg = Vm sin(ωt) si calcola la tensione di uscita raddrizzata: VDC 1 = T ˆT Vm sin(ωt)dt = 0 T 1 Vm Vm (1 − cos(ωt))|02 = T π ed il suo valore quadratico medio Vrms definito come: 2 Vrms 1 = T ˆT 0 Vm2 sin(ωt)2 dt = Vm2 Vm ⇒ Vrms = 4 2 49 Raddrizzatori Da questi valori è possibile calcolare sia la potenza in ingresso: PAC = Vrns Irms Vm V2m Vm2 = = 2 R 4R che quella in uscita: PDC = VDC IDC = Vm2 π2R ottenendo il rendimento: η= PAC π2 = PDC 4 Caso con carico R − Vx andiamo ora a vedere cosa accade se il carico è pari ad un generatore di tensione Vx . Com’è intuitivo, fintanto che non si verifica la condizione Vg > Vx (αin Figura 5.3) non si ha alcuna tensione sull’uscita. Si ottiene così: α : Vx = Vm sin(α) da cui deriva α = arcsin VDC IDC 1 = 2π 1 = 2π Vx Vm e: π−α ˆ Vm sin(ωt)dωt α π−α ˆ Vm sin(ωt) − Vx dωt R α Caso con carico induttivo R−L cosa accade se il carico invece di essere puramente resistivo ha anche una componente induttiva? Nel caso di carico R-L, durante il periodo in cui l’interruttore è attivo deve valere la relazione: Vg = L · 50 di(t) + R · i(t) dt di C.Bettini, A.Guiggiani e T.Lorenzetti Raddrizzatori Figura 5.1.: Caso con carico R − Vx quindi l’onda non è più regolare ma presenta un piccolo “strascico” che la fa finire oltre il punto π (Figura 5.2). Si possono individuare tre istanti di tempo: 1. per 0 < t < α si ha Vg > VR da cui deriva il fatto che VL > 0 e, quindi, l’induttore si carica; 2. per α < t < π si ha Vg < VR e, quindi, l’induttore si scarica; 3. per π < t < β si ha Vg < 0 ma, al contempo, l’induttore in scarica mantiene aperto il diodo e, quindi, parte della tensione negativa passa sul carico. Si ottiene così: VDC 1 = 2π π+α ˆ Vm Vm Vm sin(ωt)dωt = (1 − cos(π + α)) 6 2π π 0 dove l’uguaglianza si raggiunge nel caso in cui α = 0 e, di conseguenza, anche Vx = 0. Caso con carico L − Vx analizziamo ora il caso in cui il carico sia equivalente ad una induttanza in serie ad un generatore di tensione Vx . Deve valere: Vg = L · di(t) + Vx dt di C.Bettini, A.Guiggiani e T.Lorenzetti 51 Raddrizzatori Figura 5.2.: Caso con carico R − L da cui deriva VL = Vg − Vx dove VL è la tensione ai capi dell’induttanza. Come nel caso precedente, il diodo è forzato a condurre fino a che l’induttore non si scarica, però notiamo una seconda particolarità data dal generatore Vx (Figura 5.3): definito α l’istante in cui Vg = Vx , allora si ottiene che il valore di β (precedentemente definito nel caso di carico induttivo) è tale per cui: β =π−α con α = arcsin Vx Vm e Vm tale che Vg = Vm sin(ωt). Raddrizzatore a doppia semionda utilizzando un trasformatore con il secondario dotato di una presa a metà avvolgimento (anche detto trasformatore a presa centrale) è possibile ottenere due tensioni sfasate di 180º che possono essere singolarmente raddrizzate per mezzo di due diodi (Figura 5.4.b e Figura 5.5.b): si noti che la 52 di C.Bettini, A.Guiggiani e T.Lorenzetti Raddrizzatori Figura 5.3.: Caso con carico L − Vx tensione totale del secondario deve essere doppia rispetto a quella necessaria per il raddrizzamento ad una semionda. Raddrizzatore a ponte adottando quattro diodi disposti in configurazione a ponte di Graetz (dal nome del suo inventore, il fisico tedesco Leo Graetz, Figura 5.6), è possibile ottenere un segnale che è la somma di una semionda positiva più la semionda negativa capovolta (doppia semionda). Questa soluzione, molto usata negli alimentatori, rende molto più semplice il successivo filtraggio e livellamento della tensione fino ad ottenere una corrente continua, non richiedendo peraltro un trasformatore con doppio avvolgimento a presa centrale. Principale svantaggio di questo metodo è di avere una caduta di tensione pari a quella di due diodi in serie, quindi anche di oltre 2V : nel raddrizzare tensioni molto piccole si ha quindi una perdita e una distorsione eccessive. di C.Bettini, A.Guiggiani e T.Lorenzetti 53 Raddrizzatori Figura 5.4.: Raddrizzatore a singola (a) ed a doppia (b) semionda Analizziamo i valori di: VDC 1 = π ˆπ Vm sin(ωt)dωt = 0 2Vm π e di: 2 Vrms 1 = π ˆπ 0 Vm2 sin(ωt)2 dωt = Vm2 Vm ⇒ Vrms = √ 2 2 Una configurazione simile costituita da sei diodi permette di raddrizzare una tensione trifase impiegando tutte e tre le fasi (anche più di tre in un sistema polifase, usando un numero opportuno di diodi). 54 di C.Bettini, A.Guiggiani e T.Lorenzetti Raddrizzatori Figura 5.5.: Uscita di un raddrizzatore a singola (a) ed a doppia (b) semionda Raddrizzatore a ponte con generatore reale e carico ID analizziamo adesso il caso in cui il raddrizzatore abbia un carico ID ed il generatore utilizzato sia modellato come reale (quindi con una induttanza parassita L in serie). Consideriamo il circuito semplificato relativo al primo diodo (D1 , Figura 5.7) in cui il diodo D2 è detto “diodo di free-wheeling” in quanto consente ad ID di circolare anche quando l’interruttore D1 è in stato di OF F . Consideriamo i vari casi: • per t < 0 il diodo D1 è in OF F mentre D2 richiude il generatore su sè stesso; • per 0 < t < u il diodo D1 è in ON e VL inizia a caricarsi fino ad un valore massimo dipendente da ID1 ; di C.Bettini, A.Guiggiani e T.Lorenzetti 55 Raddrizzatori Figura 5.6.: Raddrizzatore a ponte • per u < t < π il diodo D1 è in ON mentre D2 è in OF F : la tensione Vout = Vg ed Iout = ID . La corrente del circuito satura ad iL . Si osservi Figura 5.8 per analizzare gli andamenti di tensioni e correnti nel circuito. Come si può notare, l’induttanza porta un ritardo sulla chiusura di D1 che, a sua volta, porta ad avere un valore medio in uscita minore del caso ideale. Con un po’ di calcoli si ottengono le relazioni: u = arccos 1 − ωLID Vm e: VDC 1 = 2π ˆπ Vm sin(ωt)dωt = 0 Vm ωLID − π 2π Raddrizzatore di precisione qualora il segnale da raddrizzare abbia una tensione molto bassa, la tensione di caduta del diodo non è trascurabile. Poiché la conduzione inizia solamente dopo il superamento del valore di soglia, segnali inferiori vengono del tutto soppressi. Anche oltre la soglia, la caduta di tensione è sottratta al segnale. Per ovviare a questo inconveniente, negli strumenti di misura e altri dispositivi dove sia richiesta una rettificazione precisa del segnale, si usano diodi inseriti nel circui- 56 di C.Bettini, A.Guiggiani e T.Lorenzetti Raddrizzatori Figura 5.7.: Raddrizzatore con generatore reale e carico ID to di retroazione di un amplificatore operazionale (Figura 5.9). In questo circuito l’amplificatore lavora come inseguitore, portando il valore di Vo allo stesso valore di Vi. Perché questa condizione si verifichi occorre che: Vi > Vd G dove Vd è la caduta di tensione sul diodo e G è il guadagno dell’amplificatore operazionale. Poiché solitamente G è nell’ordine delle centinaia di migliaia, la tensione di soglia è ridotta di un equivalente fattore rispetto alla tensione di caduta e quindi l’errore è dato principalmente dagli errori dell’amplificatore operazionale, in particolare causato da sbilanciamento della tensione d’ingresso, dalla velocità e dalla corrente d’ingresso del piedino invertente. Raddrizzatore trifase gli schemi proposti precedentemente sono relativi a sistemi monofase; per sistemi trifase, bisogna utilizzare un’altra tipologia di raddrizzatore il cui schema è mostrato in Figura 5.10. Il principio di funzionamento del raddrizzatore trifase è lo stesso di quello monofase; ciascun diodo conduce quando la tensione tra l’anodo e il catodo è maggiore di zero (nei diodi reali la tensione deve essere maggiore di 0.7V ) e risulta interdetto nel caso contrario. In Figura 5.11 e Figura 5.12 viene mostrato in dettaglio lo stato di ciascun diodo in funzione delle tre tensioni di ingressoe l’andamento della tensione in uscita (tra i punti P-N) in assenza del condensatore (che normalmente viene inserito in uscita per livellare la tensione e il cui dimensionamento verrà trattato in seguito). di C.Bettini, A.Guiggiani e T.Lorenzetti 57 Raddrizzatori Figura 5.8.: Uscita e andamenti di tensione e corrente in un raddrizzatore con carico ID Il valore medio della tensione in uscita al raddrizzatore trifase può essere calcolato mediante la relazione: √ 2 Add Vd = √ = √ Vpp ≈ 1.35Vpp 3 3 Dimensionamento del condensatore come si è visto precedentemente, la tensione in uscita ad un raddrizzatore presenta delle ondulazioni che devono essere necessariamente ridotte a valori accettabili mediante l’utilizzo di un condensatore collegato 58 di C.Bettini, A.Guiggiani e T.Lorenzetti Raddrizzatori Figura 5.9.: Raddrizzatore di precisione Figura 5.10.: Raddrizzatore trifase in parallelo al carico. Nell’intervallo t0 − t1 il diodo è polarizzato direttamente, il condensatore si carica con una costante di tempo molto bassa (τc = Rf C, dove Rf è la resistenza interna del diodo in conduzione). All’istante t1 il catodo del diodo si trova alla tensione Vp (tensione di picco) che è la tensione a cui si è caricato il condensatore; il diodo si polarizza inversamente e il condensatore inizia a scaricarsi con costante di tempo che dipende dal carico e dalla capacità (un esempio in Figura 5.13). Si dimostra facilmente che nel caso di raddrizzatore a singola semionda, la tensione Vrpp è data dalla relazione: Vrpp = VDC f · RC di C.Bettini, A.Guiggiani e T.Lorenzetti 59 Raddrizzatori Figura 5.11.: Stato dei diodi e uscita di un raddrizzatore trifase Figura 5.12.: Uscita di un raddrizzatore trifase mentre nel caso di raddrizzatore trifase si dimensiona il condensatore considerando: C= T I · 4 ∆V che lega la capacità stessa al periodo del segnale da raddrizzare (T ), al ripple di tensione massimo (∆V ) che si vuole ottenere ed alla corrente (I). I condensatori che si trovano in commercio hanno in genere un limite di tensione di 450V . Ciò rappresenta un grosso limite nel caso in cui il condensatore deve essere collegato a valle di un raddrizzatore trifase. In questo caso, infatti, il valore di 60 di C.Bettini, A.Guiggiani e T.Lorenzetti Raddrizzatori Figura 5.13.: Esempio di condensatore di carico tensione in uscita dal raddrizzatore è di circa 513V e il problema può essere facilmente risolto utilizzando una configurazione come quella mostrata in Figura 5.14: collegando due condensatori in serie, la tensione tra i punti P-N si ripartisce in due parti uguali; nello stesso tempo, però, anche la capacità totale dei due condensatori si dimezza. Per riottenere nuovamente la capacità totale richiesta è sufficiente collegare in parallelo ai due condensatori un’altra coppia di condensatori. In definitiva si ottiene una serie di due paralleli. Figura 5.14.: Possibile alternativa al condensatore unico di C.Bettini, A.Guiggiani e T.Lorenzetti 61 Raddrizzatori Note e parametri un ponte raddrizzatore può essere realizzato collegando opportunamente diodi semplici, ma esiste in commercio una grande varietà di dispositivi integrati pronti all’uso. La soluzione a diodi singoli è particolarmente indicata per correnti elevate, poiché facilita lo smaltimento del calore. Il raddrizzatore può ricevere la corrente alternata da un trasformatore riduttore o direttamente dalla rete elettrica. Il segnale pulsante in uscita da un raddrizzatore può essere considerato come la sovrapposizione di una componente alternata e una componente continua che ne trasla il valore medio. Per questo, per livellare la corrente continua prodotta, si pone all’uscita del raddrizzatore un circuito RC passa-basso in modo da sopprimere la componente alternata. Spesso la resistenza non è aggiunta, ma costituita dalle resistenze interne dei conduttori, dei diodi e del condensatore. Per avere un’idea di massima della tensione livellata si usa la formula: √ V0 = Vef f 2 dove Vef f è la tensione efficace in entrata al raddrizzatore. Sia nel caso di raddrizzatori ad una semionda che a due semionde, i parametri da considerare per la scelta del dispositivo sono: • tensione massima nominale: ciascun diodo deve sopportare senza guastarsi una tensione inversa pari alla tensione di picco. Comunemente sono usati diodi con tensioni di breakdown di oltre 1 kV. • corrente massima nominale: questo parametro è funzione della potenza del dispositivo usato. Si tratta della corrente sopportata con continuità supponendo che il dispositivo sia correttamente raffreddato. Sono possibili momentanee sovracorrenti, per esempio all’accensione durante la carica dei condensatori. Note sul rendimento energetico Ciascun diodo, quando è attraversato da corrente, presenta una caduta di potenziale ai suoi capi relativamente costante. Per i diodi al silicio questo valore è intorno ad 0.7 ÷ 1V . La potenza dissipata da ciascun diodo è data dalla tensione presente ai suoi capi per la corrente che lo attraversa. Poiché in un ponte raddrizzatore, durante ogni semionda, conducono due diodi, la potenza totale è pari al doppio di quella dissipata da un singolo diodo. Ad esempio, 62 di C.Bettini, A.Guiggiani e T.Lorenzetti Raddrizzatori supponendo una caduta di 1V , con una corrente efficace di 10A, ciascun diodo dissipa 10W . Poiché due diodi conducono in ogni istante, la potenza totale continua dissipata è di 20W , che vanno sottratti alla potenza in entrata per ottenere il valore di potenza erogata in uscita. È intuitivo il fatto che i dispositivi di questo tipo debbano essere raffreddati (a meno che le correnti non siano limitate a pochi ampere), generalmente per mezzo di alette metalliche ed eventualmente ventilatori. Nei più grandi apparati di raddrizzamento il raffreddamento è spesso svolto da un circuito idraulico. di C.Bettini, A.Guiggiani e T.Lorenzetti 63 6. Azionamenti elettrici Quando si parla di azionamenti elettrici si tende a suddividere l’argomento in trasformatori e motori elettrici. Per quanto riguarda motori elettrici si deve inoltre distinguere tra funzionamento asincrono e sincrono: nel primo caso i campi elettrici caratterizzanti rotore e statore ruotano in modo asincrono caratterizzato da un certo sfasamento (s), nel secondo in modo sincrono. 6.1. Cenni sui trasformatori Il trasformatore è una macchina elettrica statica (in quanto non contiene parti in movimento) e reversibile, che serve per variare i parametri della potenza elettrica apparente (tensione e corrente) in ingresso rispetto a quelli in uscita (mantenendola al contempocostante). Il rendimento di un trasformatore è molto alto e le perdite sono molto basse (nel ferro, per effetto dell’isteresi e delle correnti parassite, e nel rame, per effetto Joule). Il trasformatore è una macchina in grado di operare solo in corrente alternata sfruttando i principi dell’elettromagnetismo legati ai flussi variabili. Il trasformatore viene ampiamente usato nelle reti di trasporto dell’energia elettrica che collegano le centrali elettriche alle utenze (industriali e domestiche) ed è stato uno dei motivi principali della vittoria della corrente alternata di Tesla nella famosa “guerra delle correnti” contro Edison. Le enormi quantità di energia elettrica richieste dalla società moderna fanno sì che questa debba essere prodotta in grandi quantità presso le centrali elettriche. Un parametro utile per determinare la dimensione e la quantità di energia prodotta da una centrale è la potenza, la quale può variare dalle decine di kW (per piccole centrali 65 6.1 Cenni sui trasformatori idroelettriche o solari) alle centinaia di M W delle grandi centrali termoelettriche e nucleari. Questa energia deve essere trasportata anche per centinaia di km. La potenza elettrica è legata in maniera diretta ai parametri di tensione e corrente, secondo la formula: P = V I cos φ dove cos φ (detto fattore di potenza) è il correttivo dovuto allo sfasamento. Ciò significa che a parità di potenza aumentando la tensione V diminuisce la corrente I (e si deve mantenere più vicino possibile al valore unitario). Ciò è molto importante in quanto la corrente I genera al suo passaggio nei conduttori elettrici calore (per effetto Joule): più la corrente è alta e più calore si genera. Per ovviare a questo inconveniente bisogna aumentare la sezione dei conduttori, ma esiste un limite economico e tecnologico nel dimensionamento delle linee elettriche (legato anche al fenomeno della caduta di tensione delle linee stesse); al fine quindi di abbassare la corrente I si effettua una trasformazione aumentando la tensione V (a parità di potenza P ). Diminuendo le distanze da percorrere e la potenza da trasportare viene anche meno l’esigenza di avere tensioni alte, se a questo si associa anche l’esigenza di avere per l’uso domestico e industriale un livello di tensione compatibile con le esigenze di sicurezza: ne segue che dalla produzione alla distribuzione sono necessarie un numero adeguato di trasformazioni verso tensioni via via più basse. La macchina elettrica che si occupa di effettuare tali trasformazioni è appunto il trasformatore. A titolo di esempio, si riporta un elenco delle tensioni tipiche di esercizio degli impianti elettrici: • 230 V: tensione per usi domestici; • 400 V: tensione per uso industriale; • 8.4/20 kV (8.400 ÷ 20.000 V): tensione di esercizio delle reti elettriche di distribuzione secondaria (fino ad alcune decine di km); • 132/150/230/400 kV: tensione di esercizio delle linee elettriche di distribuzione primaria (fino a centinaia di km). 66 di C.Bettini, A.Guiggiani e T.Lorenzetti 6.1 Cenni sui trasformatori Il trasformatore più semplice è costituito da due conduttori elettrici (solenoidi) avvolti su un anello di materiale ferromagnetico detto nucleo magnetico. L’avvolgimento al quale viene fornita energia viene detto primario, mentre quello dalla quale l’energia è prelevata è detto secondario. I trasformatori sono macchine reversibili, per cui questa classificazione non corrisponde ad un avvolgimento fisico unico. Quando sul primario viene applicata una tensione elettrica alternata sinusoidale, per effetto dell’induzione magnetica si crea nel nucleo un flusso magnetico con andamento sinusoidale. Per la legge di Faraday-Neumann-Lenz, questo flusso variabile induce nel secondario una tensione sinusoidale. La tensione prodotta nel secondario è proporzionale al rapporto tra il numero di spire del primario (NP ) e quelle del secondario (NS ) secondo la relazione: NP VP = =k VS NS Definito il valore efficace E di una tensione sinusoidale di ampiezza massima EM come: EM E= √ 2 da cui deriva che, trascurando le perdite, la relazione tra tensione, numero di spire, intensità di flusso e sezione del nucleo è data dalla relazione: ω 2π E = √ N SB = √ f N SB = 4.44 f N SB 2 2 dove E è il valore efficace (RMS) della tensione indotta, f è la frequenza in Hertz, N è il numero di spire dell’avvolgimento al quale si fa riferimento, S è la sezione del nucleo (in m2 ) e B è il valore dell’induzione in Tesla. Per trasformatore ideale si assume la convenzione degli utilizzatori dal lato primario e quella dei generatori dal lato secondario. Questo è governato dalle equazioni simboliche: V1 = kV2 I1 = I2 k di C.Bettini, A.Guiggiani e T.Lorenzetti 67 6.1 Cenni sui trasformatori dove k è il "rapporto di trasformazione". Un trasformatore reale approssima quello ideale quando: la riluttanza del nucleo è nulla (quindi la permeabilità magnetica del nucleo è infinita); le perdite nel nucleo sono nulle (quindi sono nulle le perdite nel ferro per correnti parassite ed isteresi magnetica); gli avvolgimenti hanno accoppiamento perfetto (quindi non vi sono flussi dispersi); le resistenze degli avvolgimenti sono nulle (senza perdite per effetto Joule). 6.1.1. Circuito reale Analizziamo adesso il circuito reale del trasformatore: dovremo considerare delle perdite per effetto Joule negli avvolgimenti (rispettivamente R1 ed R2 ) e delle perdite di flusso (rappresentabili attraverso due induttanze L1 ed L2 ); oltre a queste saranno presenti delle perdite nel ferro (“eddy current”) rispettivamente per effetto Joule (rappresentabili come Rm ) e per isteresi (rappresentabili come Lm ). Si noti in Figura 6.1 lo schema equivalente del trasformatore reale. Figura 6.1.: Schema equivalente del trasformatore reale Notiamo adesso come è possibile trasportare le grandezze del primario sul secondario e viceversa. Ignorando il parallelo formato da Lm ed Rm , è possibile considerare Z1 = R1 + L1 e Z2 = R2 + L2 tali che: V2 = Z2 = I2 N2 V1 N 1 N1 I N2 1 N2 = N1 da cui deriva, definendo n = N2 N1 2 V1 N2 = I1 N1 2 Z1 , la relazione: Z2 = n2 Z1 68 di C.Bettini, A.Guiggiani e T.Lorenzetti 6.2 Motore asincrono che caratterizza la relazione tra componenti presenti sul primario e loro riflesso sul secondario. 6.2. Motore asincrono Il motore asincrono è un tipo di motore elettrico in corrente alternata in cui la frequenza di rotazione non è uguale o un sottomultiplo della frequenza di rete, ovvero non è "sincrono" con essa; per questo si distingue dai motori sincroni. Sono costituiti da uno statore su cui sono avvolti i solenoidi che costituiscono le espansioni polari. In un motore trifase le espansioni sono in numero multiplo di tre. All’interno dello statore si trova un rotore solidale all’albero costituito da un pacco di lamierini magneticial silicio in cui è incluso (in genere in alluminio pressofuso ) un circuito detto a gabbia di scoiattolo, poiché è costituito da una serie di barre disposte a formare un cilindro tra due anelli. È anche detto motore con rotore in cortocircuito. Quando lo statore genera un campo magnetico rotante, nella gabbia si inducono forze elettromotrici elettriche per la legge di Faraday-Neumann-Lenz, e tali forze elettromotrici sostengono correnti. La mutua interazione attraverso i relativi campi magnetici tra le correnti di rotore e quelle di statore produce una coppia risultante netta. Esiste una velocità n alla quale il motore non produce coppia, detta velocità di sincronismo. Essa non viene usualmente superata ed è legata alla frequenza f di alimentazione ed al numero di coppie polari p dalla relazione: n = 60 · f p Per esempio, un motore con tre coppie polari (6 poli totali), alimentato a 50 Hz ha una velocità angolare teorica di 1000 giri al minuto. La velocità reale è sempre minore di circa il 3 ÷ 6%: è il fenomeno dello scorrimento che consente la produzione della coppia. Dalla formula dello scorrimento è possibile esprimere la velocità di rotazione effettiva del rotore (nr ): s= ns − nr ns di C.Bettini, A.Guiggiani e T.Lorenzetti 69 6.2 Motore asincrono dove s è lo scorrimento, a significare che il rotore scorre, cioè perde giri rispetto allo statore, ns è la velocità in numero di giri al minuto del campo magnetico di statore ed nr è la velovità del rotore. Ovviamente la presenza di un carico può ridurre ulteriormente la velocità. Se s fosse uguale a 0 vorrebbe dire che il rotore sarebbe in perfetto sincronismo, cioè avrebbe la stessa velocità del campo magnetico rotante ns . Infatti se fosse nr = ns allora risulterebbe ns − nr = 0. Se, invece, lo scorrimento s fosse uguale ad 1, il rotore è fermo. Infatti, rotore fermo implica nr = 0. Quindi lo scorrimento è uguale ad 1 quando il rotore è fermo, ovvero alla partenza. Il rendimento η del motore asincrono trifase lo possiamo calcolare con la solita formula: η= Pr Pa dove Pr e Pa sono rispettivamente la potenza meccanica utilizzata sul rotore e la potenza elettrica assorbita sullo statore. La potenza sullo statore è di tipo elettrico, quindi la si può misurare con dei wattmetri. Essendo la potenza sul rotore di tipo meccanico la possiamo trasformare in potenza di tipo elettrico se si calcolano le perdite, cioè la potenza perduta Pp . Le perdite si potenza sono dovute sia al riscaldamento degli avvolgimenti di statore e di rotore, per effetto Joule, sia alle perdite nel ferro dovute ai flussi magnetici dispersi nello statore e nel rotore, e sia alle perdite dovute agli attriti meccanici e alle ventole di raffreddamento. Se indichiamo con Pp la somma di tutte le perdite, allora la potenza resa sul rotore sarà data dalla relazione: Pr = Pa − Pp Di conseguenza: η= 70 P a − Pp Pa di C.Bettini, A.Guiggiani e T.Lorenzetti 6.2 Motore asincrono Il rendimento è basso per i piccoli motori, intorno al 77%, mentre è elevato per i grandi motori (raggiunge il 97%). Questi motori sono largamente utilizzati nell’industria a causa della elevatissima affidabilità dovuta all’assenza di spazzole. Gli avvolgimenti statorici sono in genere inglobati in resine che garantiscono un’ottima protezione dall’acqua e dagli agenti atmosferici. Questi motori sono frequentemente alimentati per mezzo di inverter elettronici che possono variarne la velocità variando in modo coordinato la frequenza e la tensione di alimentazione. L’uso di inverter permette di azionare il motore anche a partire da una corrente continua, come avviene nella trazione ferroviaria. Gli avvolgimenti statorici trifase possono essere collegati a stella oppure a triangolo, permettendo di alimentare lo stesso motore con tensioni trifase di 400V e 230V . In alcuni grossi motori si preferisce avviare a stella e poi commutare a triangolo, al fine di limitare le correnti di spunto, quando non sono utilizzati gli inverter. Per quanto riguarda la struttura dei circuiti indotti, il tipo modello di rotore più semplice e robusto si realizza infilando nei canali altrettante sbarre di rame, ciascuna delle quali riempie completamente un canale. Le testate delle sbarre che sporgono dal pacco lamellare vengono direttamente collegate fra loro, da una parte e dall’altra, mediante un grosso anello di rame. Il rotore così costruito prende la forma e viene indicato coi nomi di rotore a gabbia di scoiattolo (Figura 6.2) o rotore in corto circuito.. Questi motori sono largamente utilizzati nell’industria in quanto affidabili ed economici. La caratteristica meccanica rappresenta l’andamento della coppia motrice C in funzione della velocità di rotazione del rotore nr . La caratteristica meccanica si può anche rappresentare in funzione dello scorrimento s (Figura 6.3). Ricordiamo che scorrimento s = 1 significa motore fermo, s = 0 invece che la velocità è la massima. Per ulteriori informazioni circa lo studio e la caratterizzazione di una macchina asincrona si fa riferimento al capitolo contenente l’esperienza relativa al corso in questione. di C.Bettini, A.Guiggiani e T.Lorenzetti 71 6.3 Motore sincrono Figura 6.2.: Motore a gabbia di scoiattolo in sezione 6.2.1. Motore asincrono monofase Per le piccole potenze si costruiscono dei motori asincroni monofasi, cioè quelli che utilizzano la comune tensione presente nelle abitazioni civili tra fase e neutro a 220V e frequenza di 50Hz. Un motore asincrono monofase è composto da due avvolgimenti: il primo avvolgimento, principale, è quello che funziona a regime e non è in grado di generare un campo magnetico rotante tale da far partire il motore; di conseguenza occorre un secondo avvolgimento detto di avviamento che ha lo scopo di far partire il motore sotto carico. L’avvolgimento di avviamento ha in serie un condensatore, il quale ha la funzione di sfasare di 90° la corrente dell’avvolgimento di avviamento rispetto a quella dell’avvolgimento principale. In tal modo si genera un campo magnetico rotante in grado di far partire il motore. Una volta partito il motore, l’avvolgimento di avviamento può essere staccato mediante un interruttore che commuta non appena sia raggiunta la velocità di regime (a causa della forza centrifuga). 6.3. Motore sincrono Il motore sincrono è un tipo di motore elettrico in corrente alternata in cui il periodo di rotazione è sincronizzato con la frequenza della tensione di alimentazione, solitamente trifase. 72 di C.Bettini, A.Guiggiani e T.Lorenzetti 6.3 Motore sincrono Figura 6.3.: Andamento di coppia e scorrimento rispettivamente per: A) Motore monofase B) Motore polifase a singola gabbia di scoiattolo C) Motore polifase a singola gabbia di scoiattolo a barre profonde D) Motore polifase a doppia gabbia di scoiattolo È costituito da un rotore (parte rotante solidale all’albero) su cui sono presenti diversi poli magnetici di polarità alterna creati da magneti permanenti o elettromagneti alimentati in corrente continua, e da uno statore su cui sono presenti gli avvolgimenti del circuito di alimentazione. Le espansioni polari dello statore creano un campo magnetico rotante che trascina le espansioni polari del rotore. La frequenza di rotazione è in relazione con la frequenza di alimentazione in funzione del numero di terne di espansioni polari presenti nel motore. L’avviamento di questo tipo di motore è relativamente complesso. A motore fermo, l’applicazione della tensione alternata fa si che il rotore, per effetto dell’inerzia non abbia il tempo di seguire il campo magnetico rotante, rimanendo fermo. Il motore viene quindi inizialmente portato alla velocità di rotazione per mezzo di un motore asincrono, quindi, dopo avere scollegato quest’ultimo, viene collegata la tensione di alimentazione ed inserito il carico meccanico utilizzatore. Un’altra tecnica di avviamento sfrutta la possibilità di fare funzionare temporaneamente come asincroni motori appositamente realizzati, quindi passare al modo sincrono. Se una volta di C.Bettini, A.Guiggiani e T.Lorenzetti 73 6.3 Motore sincrono a regime la rotazione viene frenata o accellerata oltre un certo limite, si innesca una serie di oscillazioni che portano il motore al blocco e possono provocare forti sovracorrenti tali da danneggiare il motore. Per questo motivo devono essere previsto un interruttore magnetotermico di protezione. A causa della limitata praticità del motore sincrono, il suo uso con alimentazione diretta dalla rete è limitato a campi di applicazione ove sia richiesta una velocità di rotazione particolarmente precisa e stabile. È invece molto usato per azionare carichi a velocità variabile ove alimentato da convertitore statico. Esistono anche piccoli motori sincroni ad avvio automatico ed alimentazione monofase utilizzati in meccanismi temporizzatori quali i timer delle lavatrici domestiche e un tempo in alcuni orologi, sfruttando la buona precisione della frequenza della rete elettrica. In alcuni casi viene anche utilizzato come sistema di rifasamento, se provvisto di avvolgimenti rotorici. Particolarità se un motore sincrono è provvisto di avvolgimenti rotorici, al variare della corrente di eccitazione che circola sul rotore, la rete elettrica di alimentazione può vedere un carico di tipo resistivo (ohmico), induttivo o capacitivo: 1. se il motore è sottoeccitato, viene visto come ohmico-induttivo; 2. se il motore è eccitato adeguatamente, viene visto come carico puramente resistivo: in questa condizione, che è quella di normale funzionamento, si ha il minimo assorbimento di corrente ed il massimo rendimento; 3. se il motore è sovraeccitato, viene visto come ohmico-capacitivo. Il motore si presenta tanto più con caratteristica ohmica tanto maggiore è la coppia resistente del carico che gli viene applicato all’albero. Se non si applica nessun carico e si alimenta il motore per ottenere la caratteristica ohmico-capacitiva, si otterrà un compensatore rotante o più propriamente un condensatore rotante, poiché la caratteristica ohmica sarà ridotta al minimo. Questo è utilizzato come sistema di rifasamento soprattutto nelle centrali di trasformazione dell’energia elettrica. La quantità di energia reattiva che può fornire il compensatore rotante è tanto maggiore quanto maggiore è la sovraeccitazione della macchina. Il compensatore sincrono, o 74 di C.Bettini, A.Guiggiani e T.Lorenzetti 6.3 Motore sincrono rifasatore rotante, è oggi perlopiù sostituito da gruppi di rifasamento composti da condensatori statici. La pulsazione di sincronismo, che è la velocità del campo magnetico rotorico e statorico, è in relazione alle coppie di poli p: ω ωs = p 2 mentre il circuito equivalente è illustrato in Figura 6.4 dove Er = kBr è la tensione indotta sul rotore (proporzionale al campo Br ) mentre Es = kBs = jXs Ia è la tensione indotta sul rotore (proporzionale al campo Bs e, di conseguenza, alla corrente di armatura Ia ). Il campo totale è dato da: Btot = Bs + Br da cui deriva la coppia sviluppata definita come: Cdev = kBs Btot sin δ dove δ è l’angolo di sfasamento tra Br e Btot (Figura 6.5). Figura 6.4.: Schema equivalente di un motore sincrono Figura 6.5.: Campi magnetici in un motore sincrono di C.Bettini, A.Guiggiani e T.Lorenzetti 75 Parte II. Affidabilità e analisi dei rischi 77 7. Definizioni e concetti di base Conformità (conformity) si tratta di “verificare la rispondenza delle caratteristiche del prodotto alle specifiche”. Si effettua tramite l’analisi di capacità di processo, tecniche statistiche usate per conoscere la percentuale di produzione che sarà conforme. Affidabilità (reliability) è la conformità nel tempo, ovvero la probabilità che un componente mantenga inalterate le prestazioni in un istante o per un determinato intervallo di tempo fissate le condizioni di impiego. Manutenibilità (maintainability) è la probabilità di un dispositivo ad essere riportato (o mantenuto) in una condizione in cui può svolgere la funzione richiesta. Disponibilità (availability) è la probabilità che il dispositivo mantenga inalterate le prestazioni nel tempo (o in un determinato intervallo) fissate le condizioni di impiego e supponendo che siano assicurati i mezzi esterni eventualmente necessari 7.1. Affidabilità La funzione di affidabilità: R(t) = e−λt è interamente caratterizzata dal tasso di guasto λ = [ore]−1 . Solitamente il tasso di guasto dei singoli componenti è contenuto in determinate banche dati o attraverso 79 7.1 Affidabilità specifiche prove di laboratorio. Sotto alcune ipotesi di indipendenza, il tasso di guasto del sistema è pari alla somma dei tassi di guasto dei componenti. È da notare una distinzione tra ambito meccanico ed ambito elettronico: nel primo caso la formula R(t) = e−λt perde di validità e si hanno pochissime banche dati (si usano molto spesso le sole prove di laboratorio); nel secondo caso, invece, un componente può essere usato su un gran numero di dispositivi diversi e si hanno quindi molti ritorni dal campo per generare banche dati consistenti. Componenti forti e componenti deboli all’uscita dal processo di produzione i componenti sono funzionanti (poichè conformi); entro poche ore i “componenti deboli” si guastano. Per non mettere in commercio tali componenti si possano fare delle prove di setacciatura (screening tests) per “uccidere” i componenti deboli subito dopo il processo produttivo: ne è un esempio il “burn-in”, tecnica che attraverso le alte temperature porta al guasto i componenti deboli. Prove di affidabilità ad esaurimento un metodo per calcolare l’affidabilità di un componente è quello di ricorrere a prove “ad esaurimento”: seguendo quanto raffigurato in figura Figura 7.1, in ogni sequenza caratterizzata da un intervallo ∆ti viene prelevato un certo numero di componenti ni0 che sarà composto da elementi sani ed elementi guasti: ni0 = ns (t) + ng (t) Figura 7.1.: Prove di affidabilità ad esaurimento Notiamo che essendo nell’ipotesi di conformità si hanno le condizioni iniziali ng (0) = 0, n0 = ns (0); si osservano più fasi (per ogni intervallo ∆ti ): 80 di C.Bettini, A.Guiggiani e T.Lorenzetti 7.2 Disponibilità 1. assestamento preliminare: inizialmente si tengono tutti i componenti alla Tamb affinchè partano tutti dalle stesse condizioni iniziali; 2. trattamento: la temperatura viene aumentata e mantenuta (prima di diminuire nuovamente) in un ambiente controllato e con specifiche maschere per un intervallo ∆ti ; 3. ri-assestamento: si riportano i campioni alle condizioni inziali di riferimento. A questo punto il procedimento viene ripetuto finchè tutti i pezzi non si sono guastati, quindi n0 (tf ) = ng (tf ) e ns (tf ) = 0. Curva a vasca Si descrive la curva a vasca in Figura 7.2: • ZONA 1: è la zona di mortalità infantile, causata da guasti per deficienza intrinseca. Non viene riscontrata nel processo di produzione e si tenta di eliminarla con tecniche di screening; ha una durata che solitamente varia tra le 100 e le 1000 ore ed è una zona che ricade nel periodo di garanzia; • ZONA 2: è la zona di vita utile, nel caso elettronico lunga solitamente centinaia di migliaia di ore. È caratterizzata dai guasti casuali (ad esempio spark di corrente) e si suppone λ(t) costante; • ZONA 3: è la zona di usura, in cui si hanno guasti per invecchiamento. Nel caso elettronico è non interessanti poichè dopo alcuni anni l’elettronica si considera superata e quindi perde di valore. Ricordiamo che nel caso meccanico la ZONA 2 è molto più stretta e, quindi, si ricorre molto più frequentemente ad azioni di manutenzione. 7.2. Disponibilità Per definire la disponibilità si definiscono tre parametri statistici: 1. MTTF = Mean Time To Failure (il componente è non riparabile); 2. MTBF = Mean Time Between Failure (il componente è riparabile); 3. MTTR = Mean Time To Restore (il componente è riparabile). di C.Bettini, A.Guiggiani e T.Lorenzetti 81 7.2 Disponibilità Figura 7.2.: Esempio di curva a vasca È da notare la differenza tra “riparare” (che indica la sostituzione del componente in seguito ad una rottura) e “restore” (che indica il riprisstino delle funzionalità di tutto il sistema, comprensivo di diagnosi del guasto, sostituzione del componente rotto e collaudo del sistema). La disponibilità asintotica del sistema si calcola come: A0 = M T BF M T BF + M T T R La funzione di disponibilità è invece definita come in Figura 7.3. Figura 7.3.: Funzione di disponibilità A(t) Parametri stocastici di affidabilità e disponibilità se un componente non è riparabile, è possibile definire il tempo medio al guasto come valore atteso del tempo al 82 di C.Bettini, A.Guiggiani e T.Lorenzetti 7.2 Disponibilità guasto: ˆ+∞ ˆ+∞ 1 M T T F = E[tg ] = t F (t)dt = R(t)dt = nel caso in cui si è in ZONA 2 = λ 0 0 La situazione cambia drasticamente nel caso di un sistema, in cui il tempo di guasto dipende da quello dei sotto-sistemi e dai componenti. Per calcolare il tempo di guasto di un sistema è necessario sapere se questo è in configurazione serie o parallelo: • configurazione serie: il sistema è funzionante solo se tutti i suoi N-blocchi sono funzionanti. Di conseguenza: Rsistema (t) = p(e1 )p(e2 |e1 ) · ... · p(eN )p(eN |e1 , e2 , ... , eN −1 ) dove p(ei ) e p(ei |ei−1 ) rappresentano rispettivamente la probabilità di guasto di ei e quella di ei supposti sani i precedenti componenti. Se supponiamo i guasti indipendenti: Rsistema (t) = N Y Ri (t) = i=1 da cui deriva che λsistema = N Y e−λi t = e−(λ1 + ... +λN )t = e−λsistema t i=1 PN i=1 λi e vale la relazione: Rsistema < min{Ri } quindi l’affidabilità del sistema è minore del minor valore di affidabilità dei singoli componenti; • configurazione parallelo: in questo caso il sistema funziona se almeno un blocco è funzionante. Si calcola la funzione di inaffidabilità (s̄): p(s̄) = p(ē1 )p(ē2 |ē1 ) · ... · p(ēN )p(ēN |ē1 , ... , ēN −1 ) = N Y i=1 fi = N Y 1 − e λi t i=1 da cui deriva che Rsistema = 1 − Fsistema > max{Ri } e, quindi, l’affidabilità del sistema è maggiore di quella del componente ad affidabilità massima; di C.Bettini, A.Guiggiani e T.Lorenzetti 83 7.2 Disponibilità • configurazione a ponte: sono particolari configurazioni in cui non è sufficiente un’analisi serie o parallelo. In questi casi si ricorre ad altri metodi di analisi (Figura 7.4): – metodo dello spazio degli eventi: si va a vedere caso per caso tutte le combinazioni (n) in cui i vari blocchi sono (o meno) guasti. Si calcola poi: Rsistema = P (n Y ) ei i=i considerando le probabilità di intersezione nel caso in cui i guasti non siano indipendenti. – metodo delle unioni: si individua il minimo insieme di percorsi, ovvero l’insieme in cui un percorso di ordine n non contiene percorsi di ordine inferiore. Nell’esempio abbiamo i percorsi D-C, A-E ed A-B-C: si calcolano le varie p(ui ) e si ottiene: Rsistema = P (u1 ) + P (u2 ) + ... + P (un ) – metodo dei tagli: si cerca un insieme minimo dei n tagli ti , dove il “taglio” è definito come insieme di elementi che se guasti causano l’avaria del sistema. In esempio si identificano i tagli A-D, A-C, C-E, B-C-E. Si identifica così: Rsistema = P (t1 ) + ... + P (tn ) – metodo della probabilità condizionata: si cerca la probabilità dei casi in cui il sistema S sia funzionante dato il sottosistema A,B,C,D od E funzionante. In questo metodo è possibile così ricavare la probabilità inversa in cui, dato il sottosistema non funzionante, il sistema non è in avaria. 84 di C.Bettini, A.Guiggiani e T.Lorenzetti 7.2 Disponibilità Figura 7.4.: Metodi di analisi per configurazioni a ponte Valutazione della disponibilità la valutazione della disponibilità può essere effettuata con tecniche quantitative (che consentono di determiunare un modello matematico A(t) in modo da valutare la disponibilità asintotica A0 ) o con tecniche qualitative, che sfruttano determinate tabelle quali FMEA, FMECA ed FTA. • FMEA: è il Failure Modes & Effects Analysis: è una analisi di tipo “bottomup”, ovvero che parte dallo studio dei sottosistemi per ricavare poi la disponibilità del sistema. Si distinguono varie fasi: – scelta del livello di analisi; – scelta del componente da analizzare; – identificazione del modo di guasto (ovvero dell’evidenza del guasto); – determinazione del meccanismo di guasto (processo che porta al guasto); – determinazione della causa del guasto (ad esempio vibrazioni o spike di corrente); – effetti del guasto, che a loro volta si distinguono in locali (come il guasto influenza il sottosistema in cui avviene) o globali (come il guasto si riflette sul funzionamento del sistema). • L’analisi FMEA ha una seconda variazione, detta FMECA (C=Criticality) in cui si prendono in considerazione altri aspetti: – occurrance: indicatore quantitativo che caratterizza la frequenza del guasto e varia da 1 a 10; di C.Bettini, A.Guiggiani e T.Lorenzetti 85 7.2 Disponibilità – severity: indice della pericolosità del guasto, va da un minimo di 1 (nessun effetto) a 10 (hazardous); – detection: indice della possibilità di diagnosticare il modo di guasto, è un indice di prevenzione; – RPN: è il Risk Priority Number ed è il prodotto dei tre indici precedenti. Bisogna tener conto dei tre valori in quanto RPN uguali possono essere generati da diverse combinazioni; – Design Review: è un elenco di azioni correttive volte ad abbassare il valore dell’RPN; – costi AP/AC: è una valutazione dei costi delle azioni preventive e correttive (non è un passo obbligatorio dell’FMECA). • FTA: è il Fault Tree Analysis: è un’analisi di tipo “top-down” previsionale, utilizzata per stimare la probabilità di accadimento di un evento (top event) basandosi sull’individuazione delle cause del guasto e considerando anche le possibili relazioni di legame. Ha due fini: – individuare tutti i possibili guasti che causano un’avaria del sistema; – analizzare le dipendenze tra i guasti di alto e di basso livello Nello studio FTA si disegna un diagramma caratterizzato da particolari simboli interconnessi: • il top event è rappresentato da un cerchio; • con un rombo si rappresentano gli eventi non sviluppati per conseguenze non rilevanti o informazioni instufficienti ma che dipende comunque statisticamente da altri eventi; • con un rettangolo si rappresentano gli eventi intermedi che si verificano come conseguenza di una o più cause combinate attraverso porte logiche (AND, OR, ...). 86 di C.Bettini, A.Guiggiani e T.Lorenzetti 7.3 Prove 7.3. Prove La qualità di un prodotto è quantificabile attraverso la sua conformità e la sua affidabilità. Le prove di qualifica permettono di sottoporre le apparecchiature a condizioni nominali limite simulando in laboratorio le condizioni operative dei dispositivi: sono prove che devono essere uniformi e riproducibili, spesso basate su norme specifiche. Le fasi di una prova sono: 1. assestamento preliminare; 2. controlli e misure iniziali; 3. trattamento; 4. riassestamento; 5. controlli e misure finali. Queste prove si suddividono inoltre in base alla loro tipologia: si parta di prove combinate (si applicano due o più sollecitazioni contemporanee), composite (si applicano due o più sollecitazioni in rapida successione) od in sequenza (si applicano due o più sollecitazioni consecutive). Spesso vengono effettuate delle “prove di sviluppo” (quindi su pochi campioni) per ottenere più informazioni possibile sul campione prima che questo venga danneggiato. 7.3.1. Prove di affidabilità le prove di affidabilità possono essere fatte in laboratorio (altamente riproducibili ma costose) o in campo (che portano solitamente a risultati più realistici ma sono poco controllabili e riproducibili). Le prove di affidabilità si distinguono in: 1. prove di conformità dell’affidabilità, in cui si verifica l’affidabilità di un dispositivo in un intervallo di tempo ben preciso; 2. prove di vita: sono prove di tipo distruttivo che a loro volta si suddividono in prove accelerate o di lunga durata. di C.Bettini, A.Guiggiani e T.Lorenzetti 87 8. Analisi dei rischi È un’analisi qualitativa eseguita durante le prime fasi della progettazione. Aiuta ad identificare i rischi potenziali (e le relative conseguenze), a scegliere le misure per trattare i rischi e garantire che i sistema sia sicuro, permette l’implementazione di modifiche durante il progetto (quindi più facili e meno costose) e di ridurre in fase progettuale i rischi e gli imprevisti. Notiamo la differenza tra pericolo (che è una potenziale fonte di danno, ovvero un insieme di condizioni che possono portare ad un incidente) ed un rischio, definito come la combinazione della probabilità del verificarsi di un danno e la sua gravità. Le fasi dell’analisi sono: 1. identificazione dei rischi; 2. determinazione di causa-effetto; 3. determinazione della probabilità che un rischio provochi un incidente; 4. scelta dei requisiti procedurali e di progetto per eliminare/controllare i rischi. Dal punto di vista matematico il richio è il prodotto tra frequenza e criticità, ed è composto da una componente deterministica (Rid ) ed una aleatoria (RN ON id ). Ricordiamo che la stima del rischio non può essere univoca, quindi è necessario adottare degli standard e delle convenzioni riconosciute e valutare l’accettabilità del rischio (circa 10−6 nei paesi occidentali) e la sua percezione. I rischi possono essere stimati con tre diversi parametri che non sono tra loro equivalenti (bisogna scegliere di caso in caso quale parametro adottare): 1. indici di rischio: sono singoli valori numerici presentati in forma tabellare e quindi non molto rappresentativi in quanto altamente generici. Sono un 89 Analisi dei rischi esempio il FAR (Fatal Accidental Rate, stima degli eventi letali per 108 ore di esposizione), l’ARD (Average Rate of Death, numero medio di decessi per unità di tempo) e l’ESCI (Equivalent Social Cost Index, costo collettivo equivalente che non tiene conto degli eventi catastrofici); 2. rischio puntuale: è pari al valore di frequenza con cui, in un determinato punto geografico, si può verificare il danno di riferimento. Il rischio puntuale è rappresentato con curve iso-rischio o con determinati profili di rischio e può essere calcolato localmente (riferito ad un bersaglio sempre presente in un dato punto geografico e senza possibilità di fuga o protezione) od individualmente (riferito ad un bersaglio reale, presente in modo discontinuo e con possibilità di protezione). In entrambi i casi si basa sull’ipotesi che i contributi di tutti gli esiti incidentali siano additivi e tiene conto di: a) natura del danno; b) probabilità del danno; c) periodo di tempo in cui il danno può manifestarsi. 3. rischio diffuso: è una misura globale del rischio a cui una collettività può essere soggetta in una determinata area. Simile al rischio puntuale, tiene conto anche di: a) tipologia della popolazione; b) presenza di bersagli particolarmente vulnerabili; c) elementi mitigatori; d) fluttuazione di presenza (giornaliera, settimanale, stagionale, ...). Il rischio diffuso si calcola principalmente in due modi: 1. curve F-N: sono curve che registrano la frequenza F di incidenti che possono determinare un danno uguale o superiore ad un dato valore Nmax . Sono diagrammi in scala logaritmica in cui è riportata la frequenza sull’asse y ed il numero di vittime su x; 2. istogrammi I-N: complementari alle curve F-N, questi istogrammi rappresentano il numero N di persone soggette ad un certo livello di rischio individuale 90 di C.Bettini, A.Guiggiani e T.Lorenzetti 8.1 Tecniche di analisi dei rischi I. In pratica sono diagrammi nei quali la popolazione è suddivisa in classi caratterizzate da diversi valori di rischio. In seguito a queste definizioni è possibile definire l’accettabilità del rischio come “il costo che si è disposti a sostenere a fronte di alcuni eventi non determinati”. Fondamentalmente è funzione di danno probabile e vantaggi conseguiti. Per valutare l’accettabilità si analizzano: • rischi-costi-benefici: si analizzano i costi ed i benefici considerando un termine correttivo negativo che rappresenta il rischio; • preferenze rivelate: si esegue una indagine di opinione e le preferenze vengono misurate, standardizzate ed usate come indice di accettabilità; • preferenze espresse: si usano indici di accettabilità già usati in passato; • standard naturali: ci si basa su studi geologici, biologici, ..., tratti da studi scientifici ed organi internazionali. È fondamentale durante l’analisi dell’accettabilità del rischio considerare gli aspetti di: • rischio oggettivo: basato su ipotesi certe e valutazioni probabilistiche, è un concetto tecnico; • rischio percepito: basato su ipotesi non certe e valutazioni personali, è un concetto profano. I fattori che influenzano accettabilità e percezione del rischio sono: numero di persone coinvolte; confronto con la mortalità naturale; obbligatorietà del rischio; irreversibilità delle conseguenze; beneficio reale (o supposto) che deriva dal rischio; livello culturale e sociale delle persone coinvolte. 8.1. Tecniche di analisi dei rischi • Preliminary Hazard Analysis: analisi preliminare (quindi basata su poche informazioni) a seguito della quale si opera una stima delle conseguenze dei top-events identificati riportando pericoli, cause, effetti e categoria; di C.Bettini, A.Guiggiani e T.Lorenzetti 91 8.1 Tecniche di analisi dei rischi • Analisi “what-if”: tecnica in brainstorming (quindi altamente dinamica) a cui si affianca il metodo delle check-list; • HazOp (Hazard and Operability Studies): si identificano le modalità di guasto e si valutano secondo una logica deviazione-causa-conseguenza-rimedio. Il sistema viene diviso in blocchi (nodi) che vengono legati con parole guida tipo “di più”, “di meno”, “invece di”, ... È una analisi che identifica le deviazioni in modo chiaro e schematico, migliore di una FMECA ma più complessa; • ETA (Event Tree Analysis): è un’analisi induttiva simile ad FTA anche se gli eventi sono analizzati attraverso una analisi binaria. Ad esempio il top-event “presenza di un incendio” si suddivide in sotto-rami ognuno caratterizzato da una sua probabilità: – fuoco riconosciuto? SI ⇒L’allarme funziona? SI ⇒ Gli innaffiatori funzionano? SI ⇒ danno limitato NO ⇒ danno esteso NO ⇒ Gli innaffiatori funzionano? SI ⇒danno limitato NO ⇒possibile fatalità NO ⇒ danno esteso • FMEA\FMECA: già analizzate in precedenza, sono tecniche utili nello studio di sistemi meccanici, elettronici o produttivi. Le analisi FMEA/FMECA sono da utilizzarsi a inizio progettazione, durante modifiche e aggiornamenti, se il budget è limitato oppure a sostegno ad una analisi FTA; • FTA: già analizzata in precedenza, si usa tanto in fase di progetto quanto in fase avanzata, nel caso in cui si verifica un guasto o per valutare determinate modifiche. Per quanto riguarda il “dove” viene applicata, l’analisi FTA si usa su impianti con guasti molto critici o caratterizzati da eventi ricorrenti. 92 di C.Bettini, A.Guiggiani e T.Lorenzetti 8.1 Tecniche di analisi dei rischi Confronto FMECA ed FTA Fondamentalmente FMECA è un’analisi adatta per sistemi: • nuovi, per i quali gli effetti dei guasti non sono noti; • senza ridondanze; • in cui non si considerano i guasti multipli. Di contro, l’analisi FTA: • necessita che gli eventi siano noti e ben definiti; • è applicabile a sistemi ridondanti; • considera i guasti multipli. In generale, quindi, la soluzione ottima è quella di applicare l’analisi FTA ad i risultati di un’analisi FMECA. In questo modo si riesce a caratterizzare in maniera completa il sistema (ed i guasti) e ad analizzare gli effetti ricorrenti (anche se derivano da modi di guasto diversi). di C.Bettini, A.Guiggiani e T.Lorenzetti 93 9. Modelli e meccanismi di guasto I meccanismi di guasto sono i processi fisico-chimici che portano al guasto del componente. È da notare la distinzione tra guasto (inteso come cessazione dell’attitudine di un dispositivo di svolgere la funzione per cui è stato progettato) e guasto reale (che è un concetto legato a “come è fatto il componente”, ovvero una correlazione tra sollecitazione, sforzo e guasto). I modelli che descrivono i meccanismi di guasto si suddividono in: 1. modello sollecitazione-resistenza: è un modello senza memoria tipico dei materiali come, ad esempio, una barra in acciaio od una molla; 2. modello degradazione-tollerabilità: è un modello con memoria tipico dei componenti elettrici e legato alla curva a vasca. È un sempio di modello di questo tipo quello di Arrhenius; 3. modello richiesta-risposta: in un modello del genere non è detto che guasto e failure avvengano insieme; una volta che si verifica un guasto ci se ne rende conto (failure) solo quando viene richiesto il servizio relativo al guasto; 4. modello tolleranza-requisito: è un modello in cui il guasto non avviene finchè si resta all’interno di una zona di tolleranza definita come spazio tra un limite superiore di specifica (LSS) ed un limite inferiore di specifica (LIS). Un’altra suddivisione dei modelli di guasto è in: • Guasti per sovraccarico: – duttilità-fragilità (Figura 9.1): legata al comportamento elasto-plastico dei materiali, su un grafico tensione/deformazione si cerca il punto in cui 95 Modelli e meccanismi di guasto il materiale arriva alla rottura analizzando un “provino”. Si definisce: σ = tensione = p̄ A0 con p̄ il carico applicato ed A0 l’area della sezione del provino e: = def ormazione = ∆l l0 con ∆l l’allungamento ed l0 la lunghezza del provino. Notiamo nel grafico che se siamo nel tratto OA togliendo il carico il materiale torna alle condizioni iniziali (non ha memoria) e la pendenza della curva dipende dal modulo di Young. Poi, nel tratto AB il materiale subisce una deformazione permanente (è una curva elasto-plastica con memoria) che va poi a peggiorare nel tratto BC fino alla rottura del provino (in pratica dopo il tratto AB il materiale è da buttare). Si considera un materiale “fragile” se questo non sopporta le sollecitazioni, ovvero se i tratti OA e BC sono di lunghezza nulla. – carico di punta: è la tensione meccanica che si genera nei punti in cui i coefficienti fisici sono diversi (ad esempio dove si ha una discontinuità del materiale). Ad esempio all’interno di un chip nei punti angolosi si ha un possibile distacco della metallizzazione che può portare ad un guasto; – scariche elettriche distruttive: possono avere cause sia “umane” che “fisiche” (ad esempio in ambienti molto secchi), entrambe portano ad un guasto catastrofico; – sovratemperatura: può essere causata da agenti esterni (ad esempio spike di corrente) o fenomeni di usura; – sfaldamenti inferfacciali tra piani cristallografici: un piano cristallografico ha normalmente una struttura molto compatta, ma tra piano e piano la struttura è più debole, quindi è possibile che vi siano degli sfaldamenti interfacciali. I materiali più soggetti a queste problematiche sono la grafite, la mica (usata molto come isolante) e tutte le deformazioni in cui si perde l’elasticità. 96 di C.Bettini, A.Guiggiani e T.Lorenzetti Modelli e meccanismi di guasto Figura 9.1.: Esempio di grafico tensione/deformazione • Guasti per usura: – corrosione: – interdiffusione: si definisce diffusione l’attitudine di una particella atomica di “migrare in uno stato di massa di un’altra sostanza”. Mentre l’interdiffusione avviene in entrambe le direzioni, la diffusione è unilaterale. Nei solidi la diffusione dipende principalmente dalla struttura molecolare, quindi dalla temnperatura e dall’energia di attivazione. I meccanismi di diffusione si classificano in: ∗ meccanismo reticolare: avviene se i due materiali hanno caratteristiche reticolari diverse. Fissata l’energia di attivazione, queste modalità dipendono fondamentalmente dalla temperatura e si suddividono in: · meccanismo per vacanza: l’atomo A migra in B ed occupa una posizione vacante; · meccanismo interstiziale: la struttura di B è allentata ed alcuni atomi di A possono prendervi parte pur non essendoci vacanze; · meccanismo di scambio: alcuni atomi di A si scambiano con gli atomi di B. ∗ diffusione per linee di dislocazione: in alcuni materiali ci sono piani cristallografici compatti che danno origini a superfici di disloca- di C.Bettini, A.Guiggiani e T.Lorenzetti 97 Modelli e meccanismi di guasto zione. Queste sono zone in cui vi è minore compattezza e possono nascere anche a causa di impurità dei materiali. La migrazione è principalmente funzione della temperatura e della densità di dislocazioni; ∗ diffusione per contorni granulari: avviene nel caso in cui B sia un materiale a grani grossi con legami deboli; si possono così creare dei percorsi interessati dal passaggio di atomi di A in B (dipende dall’energia di attivazione e dalla temperatura); – propagazione di una spaccatura per fatica: nel caso in cui, ad esempio, il materiale presenti una piccola crepa che si propaga a causa di fattori climatici; – logorio: è un problema molto sentito, ad esempio, nei connettori. Se infatti il punto di connessione si usura, la superficie di contatto diminuisce e al contempo aumenta l’intensità di corrente (e la temperatura del materiale); in questo modo aumenta nuovamente l’usura e il pezzo alla fine si guasta. • In base alle sollecitazioni applicate: – guasti meccanici; – guasti termici; – guasti elettrici; – guasti per radiazioni; – guasti chimici. In generale i guasti per sollecitazione esterna rappresentano il 50 − 60% dei guasti totali e si suddividono in: • EOS (Electrical Over Stress): sono guasti per uno spike di corrente od una sovratensione, avvengono principalmente durante l’accensione o lo spegnimento; • ESD (Electrical Static Discharge): derivano dal fenomeno della tribo-elettricità, più accentuato in ambiente secco, che nasce quando sostanze diverse entrano 98 di C.Bettini, A.Guiggiani e T.Lorenzetti Modelli e meccanismi di guasto in conttatto tra loro. I modelli ESD sono classificati da classe 0 (0-250V) a classe 3B (>8000V). Si riportano due modelli ESD (Figura 9.2): – modello human-body: simula il comportamento di un operatore che sta lavorando ad un componente; si hanno due circuiti, quello di carica e quello di scarica; – modello charged-device-model: è usato per simulare un circuito integrato carico che tocca un punto a bassa impedenza verso massa. Figura 9.2.: Modelli ESD human body (a) e charged devide model (b) Per quanto riguarda le protezioni contro i guasti ESD, queste si classificano in: 1. protezioni attive, ad esempio l’inserimento di diodi nel circuito per proteggerne alcune parti, sono protezioni che riguardano il costruttuore; 2. protezioni passive, ad esempio tavoli o braccialetti elettrostatici, sono protezioni che riguardano l’operatore. Durante l’analisi di guasti ESD, inoltre, è necessario tener conto di trasparenza, robustezza, prontezza (data la rapidità degli eventi ESD) e compatibilità. di C.Bettini, A.Guiggiani e T.Lorenzetti 99 10. Safety e Security Entrambe questi termini significano “salvaguardia dell’incolumità delle persone” ma sono caratterizzati da diverse cause: • safety è la salvaguardia da eventi indipendenti da precise volontà, quindi accidentali; • security è la salvaguardia da attacchi o danni contro la persona (o i beni) perpetrati volontariamente da individui o gruppi. Analizzando la safety si nota che non è possibile avere un sistema fisico con errori nulli, che non esistono operazioni a rischio zero e che l’errore umano non può essere eliminato completamente. Dal 1970 (il boom economico) ad oggi sono state definite alcune componenti che aiutano a descrivere la safety: • functional safety: è la parte di safety in cui sistemi o apparati operano in risposta ad input (ad esempio se aumenta la pressione un sensore indica l’evento e spegne il sistema); • safety function: sono funzioni adottate dal sistema strumentale di sicurezza (SIS) per mantenere o riportare il sistema (o il processo) in sicurezza dopo un evento pericoloso; • safety related systems: sono i sistemi che implementano delle safety functions. Si definisce così un indice di Safety Integrity dato dalla sottrazione di rischio esistente e rischio tollerabile. Oltre a questo sono stati opportunamente definiti dei Safety Integrity Level (SIL), rappresentazioni statistiche dell’affidabilità di un sistema quando viene effettuata una richiesta di processo. Solitamente sono previsti solo tre SIL (tabella Tabella 10.1), però ne è stato recentemente introdotto un quarto (SIL4) che caratterizza i sistemi super-affidabili. 101 Safety e Security SIL 4 3 2 1 Malfunzionamenti pericolosi / ora Malfunzionamenti su richiesta (PFD) 10−8 − 10−9 10−4 − 10−5 10−7 − 10−8 10−3 − 10−4 10−6 − 10−7 10−2 − 10−3 −5 −6 10 − 10 10−1 − 10−2 Tabella 10.1.: SIL Solitamente si richiede un SIL1 nel caso di produzione e proprietà primaria, un SIL2 nel caso di produzione e proprietà secondaria (con possibili lesioni fisiche), un SIL3 nel caso di protezione di dipendenti/comunità ed un SIL4 nel caso di protezione da eventi catastrofici. Si può così definire la struttura dei Safety Instrumented System suddivisi in tre sottosistemi: 1. rilevatori, ovvero i sensori che generano le variabili da controllare; 2. l’unità di controllo (solitamente un PLC), che gestisce le connessioni ed elabora; 3. i dispositivi di blocco, che agiscono nel caso in cui si manifesti un problema di sicurezza. I sistemi SIS tengono in considerazione i guasti sicuri (ovvero i casi in cui il sistema va in blocco in modo non necessario) ed i guasti pericolosi (se il sistema è in presenza di una situazione di pericolo e non riesce a portarsi in uno stato di sicurezza) e sono quantificabili attraverso due indici: • Diagnostic Coverage: è il rapporto tra guasti casuali hardware dannosi ma rivelabili (attraverso l’auto-diagostica) e guasti casuali hardware totali; • Safe Failure Fraction: è il rapporto tra la somma di guasti hardware sicuri e dannosi ma rivelabili ed il totale dei guasti hardware. 102 di C.Bettini, A.Guiggiani e T.Lorenzetti Parte III. Sistemi Digitali 103 11. Conversione Analogico-Digitale (ADC) 11.1. Campionamento e quantizzazione Nella teoria dei segnali il campionamento è una tecnica che consiste nel convertire un segnale continuo nel tempo in un segnale discreto, valutandone l’ampiezza a intervalli di tempo regolari. In questo modo, a seguito di una successiva operazione di quantizzazione e di conversione, è possibile ottenere una stringa digitale (discreta nel tempo e nell’ampiezza) che approssimi quella continua originaria. In parole povere il campionamento consiste nell’andare a misurare il valore del segnale analogico in diversi istanti di tempo. Il tempo Tc che intercorre tra una valutazione e l’altra si chiama periodo di campionamento. La frequenza di campionamento fc = Tc−1 è invece il reciproco del periodo di campionamento. Definita fBmax la massima frequenza dello spettro del segnale da campionare, se viene rispettata la condizione di Shannon: fc > 2fBmax allora è possibile ricostruire, con l’utilizzo di apposite funzioni interpolatrici, il segnale analogico senza perderne alcuna informazione; generalmente per una buona e fedele ricostruzione del livello analogico è richiesta una frequenza di campionamento che sia dalle 5 alle 10 volte superiore alla frequenza massima contenuta nel segnale campionato. Qualora invece non venga rispettata tale condizione, si riscontra un effetto conosciuto con il nome di Aliasing, che comporta una distorsione del segnale 105 11.1 Campionamento e quantizzazione analogico ricostruito rispetto a quello originale campionato. Definendo la frequenza di Nyquist come fn = f2c , la ricostruzione di un segnale analogico costituito da componenti in frequenza superiori alla frequenza di Nyquist porta ad ottenere un segnale dove tali componenti hanno una frequenza detta "specchiata" rispetto la frequenza di Nyquist (cioè simmetrica rispetto alla frequenza di Nyquist a quella del segnale analogico originale). Ad esempio se il livello analogico è una sinusoide di frequenza fB = 12Hz e la frequenza di campionamento è fc = 20Hz, allora si ottiene fr = fn −(fs −fn ) = 8Hz, con fr la frequenza della sinusoide ricostruita. Per evitare il fenomeno dell’aliasing è necessario: 1. adottare una frequenza di campionamento superiore (se non si vogliono perdere le informazioni contenute nelle componenti ad alta frequenza del segnale analogico acquisito); 2. adottare un filtro anti-aliasing (passa-basso) in modo tale da eliminare le frequenze contenute nel segnale analogico superiori alla frequenza di Nyquist del campionatore. Quando si parla di sistemi di conversione Analogico-Digitale (ADC, figura Figura 11.1) si suddivide così il problema in: Figura 11.1.: Convertitore Analogico/Digitale (ADC) 1. campionamento, che può essere di tipo: a) over-sampling: se fc fB , è il campionamento più usato poichè molto comodo se a valle aggiungiamo un filtro (ad esempio un COMB+FIR) ed un decimatore che riporta la frequenza a quella desiderata; b) campionamento alla fsh di Shannon fc = 2fB , ideale; 106 di C.Bettini, A.Guiggiani e T.Lorenzetti 11.1 Campionamento e quantizzazione c) under-sampling: se fc fB , è usato nei segnali a “banda stretta” (se la portante f0 è molto maggiore della banda del segnale, quindi f0 fB ). 2. conversione a) quantizzazione: durante la quantizzazione il segnale può essere approssimato per: i. troncamento: comporta un errore quadratico medio pari a Q err = √ 3 con Q =passo di quantizzazione; ii. arrotondamento: comporta un errore quadratico medio pari a Q err = √ 12 b) digitalizzazione Se analizziamo la conversione notiamo che il tempo di conversione Tc è tale che: Tc 6 ∆t = 1 π2N fc con N =numero di bit del segnale ed fc frequenza di campionamento. Sample & Hold il sample & hold è uno strumento utilizzato per campionare il segnale e mantenerlo costante durante il cosìddetto conversion-time. Per fare ciò, questo apparecchio può lavorare: • in circuito aperto: con questa configurazione (figura Figura 11.2.a) si massimizza la velocità di campionamento a discapito dell’impossibilità di riconoscere eventuali offset o disturbi del segnale in ingressi. Ogni componente ha una capacità di accoppiamento che si fa sentire in modo maggiore all’aumentare della frequenza; di C.Bettini, A.Guiggiani e T.Lorenzetti 107 11.1 Campionamento e quantizzazione • in circuito chiuso: con questa configurazione (figura Figura 11.2.b) la retroazione diminuisce da un lato la banda passante (costringendoci ad adottare frequenze di campionamento minori) massimizzando però l’accuratezza. Figura 11.2.: Schema equivalente di un Sample & Hold I parametri che caratterizzano un S&H si suddividono in: • parametri statici: – durante l’azione “sample”: sono la banda passante (analogue bandwith) e la pendenza massima (max-slew rate); – durante l’azione di “hold”: sono il feedthrough (il “trasudamento” dell’ingresso sul segnale di uscita, che comporta delle oscillazioni sull’uscita) ed il droop-rate (sono delle perdite di corrente causate dai componenti attivi, principalmente gli operazionali, che fanno sì che l’uscita non sia % ). costante ma tenda a decrescere con un certo valore misurato in ns • parametri dinamici: – durante il passaggio da “hold” a “sample”: si deve attendere un tempo detto “acquisition-time”; – durante il passaggio da “sample” a “hold” si definiscono tre parametri: 108 di C.Bettini, A.Guiggiani e T.Lorenzetti 11.1 Campionamento e quantizzazione ∗ aperture delay: ritardo nell’apertura degli interruttori; ∗ aperture jitter: rumore di apertura proporzionale con la frequenza; nel caso di segnali a banda stretta questa quantità implica un ritardo descrivibile dalla formula tAj = 1 2π N f Bmax con fBmax la frequenza massima del segnale (pari a f0 + fB ); ∗ settling time: tempo necessario affinchè il segnale si assesti. Limiti di un Sample & Hold analizziamo adesso i limiti di un Sample & Hold, definendo la frequenza massima di campionamento come: fcM AX = t Ac 1 + tAp + tset + th definito il tempo di acquisizione: taq = tAc + tAp + tset che comprende l’apertura dello switch e la carica del condensatore, e il tempo di hold: tconvmin < th < thmin La scelta del S&H implica un vincolo sulla frequenza del segnale massima raggiungibile. Questa è pari a: fc SRmax 1 = min , BSH , , 2 2πA π2N ∆t ( fBmax ) con BSH pari alla banda del Sample & Hold, SRmax lo slew-rate massimo ed A l’ampiezza del segnale. Analizziamo nel dettaglio le grandezze (Figura 11.3): di C.Bettini, A.Guiggiani e T.Lorenzetti 109 11.1 Campionamento e quantizzazione • acquisition time: è il ritardo tra l’attivazione del comando di sample e l’istante in cui la tensione del condensatore è pari a quella di ingresso; • aperture time: è l’intervallo di tempo tra l’istante di attivazione dello stato di hold e quello in cui il condensatore risulta scollegato dall’ingresso, quindi il valore catturato corrisponde al valore finale dell’aperture time; • aperture jitter: intervallo di variabilità del tempo di apertura, stabilisce la massima frequenza di campionamento consentita al Sample & Hold; • setting time: è l’intervallo di tempo tra l’attivazione dell’hold e l’istante in cui l’uscita si stabilizza sul valore effettivo. Figura 11.3.: Caratteristiche e parametri di un Sample & Hold Notiamo infine che aumentando la capacità di hold si va a limitare (da un lato) il massimo slew-rate della fase di sample garantendo però una costante di tempo più alta e, di conseguenza, un minor decadimento dell’uscita nella fase di hold. Analizziamo adesso una serie di strutture di conversione partendo da sistemi semplici fino ad arrivare a casi più complessi. 110 di C.Bettini, A.Guiggiani e T.Lorenzetti 11.2 Convertitori 11.2. Convertitori Convertitori a rampa sono convertitori (Figura 11.4) abbastanza lenti in cui la relazione fra frequenza di clock (fck ) e frequenza di conversione massima (fc ) dipende strettamente dal numero di bit (N ): fc = fck 2N in pratica il convertitore a rampa va ad incrementare il valore della conversione fintanto che questo è minore dell’input analogico. Quando il comparatore ci informa che siamo arrivati al giusto livello il counter viene bloccato e la conversione è completa. Figura 11.4.: Convertitore a rampa Convertitori ad approssimazioni successive (SAR) sono convertitori (Figura 11.5) in cui la frequenza di conversione massima è pari a: fc = fck N quindi più rapidi dei convertitori a rampa. Questi convertitori approssimano il valore di conversione andando a variare i bit (dal meno significativo al più significativo) e confrontando il segnale in uscita (digitale, convertito in analogico con un DAC) con di C.Bettini, A.Guiggiani e T.Lorenzetti 111 11.2 Convertitori quello di riferimento. In questo modo con una sola “passata” (quindi in N-passi, dove N è il numero di bit a disposizione) si ottiene una conversione del valore scelto. Figura 11.5.: Convertitore ad approssimazioni successive (SAR) Convertitore Sigma-Delta Σ∆ sono convertitori più complessi e al contempo molto utilizzati nelle applicazioni reali per il costo molto contenuto. Di base i convertitori sigma-delta assicurano alta risoluzione e basso costo a scapito di basse frequenze d’ingresso. Il convertitore sigma-delta (Figura 11.6) è composto da un modulatore (Figura 11.7 e Figura 11.8) e da un filtro passa-basso (solitamente un FIR, Figura 11.9) a cui segue un decimatore (dato che la conversione viene fatta in over-sampling il decimatore ci aiuta a tornare a frequenze minori (più semplici da elaborare nel dominio digitale). Il modulatore sigma-delta lavora sempre sovra campionando il segnale: in questo modo, però, il rumore è pesato in frequenza: si parla di fenomeno di noise-shaping in quanto a basse frequenze predomina il segnale mentre a frequenze più alte predomina il rumore di quantizzazione (per questo viene filtrato il tutto). Quindi il convertitore sigma-delta è tale per cui: • in uscita l’alternanza di uni e zeri dipende da Vin ; • la frequenza (fck ), è pari alla frequenza di campionamento; 112 di C.Bettini, A.Guiggiani e T.Lorenzetti 11.2 Convertitori • la retroazione forza il valor medio del segnale logico in uscita ad eguagliare la tensione Vin ; • il rumore di quantizzazione è filtrato alle basse frequenze. A seguito del filtraggio si va a decimare il segnale in modo tale da ricondursi ad una frequenza minore e compatibile con gli apparecchi a valle del convertitore. La giusta sequenza è quindi: oversampling - filtraggio - decimazione. Figura 11.6.: Convertitore sigma-delta (Σ∆) Figura 11.7.: Modulatore sigma-delta (Σ∆) Andiamo adesso ad analizzare i convertitori ad alta velocità. Convertitori flash sono convertitori (Figura 11.10) in cui il segnale viene codificato in un singolo colpo di clock: per fare ciò vengono sfruttati un numero di comparatori pari alla risoluzione massima (nell’esempio 2N = 28 = 256) che tramite delle partizioni resistive vanno a determinare i bit del segnale. In questo modo è possibile arrivare a frequenze di conversione nell’ordine dei 5M Sps (5 MegaSamples per second) ma si è vincolati a un numero di bit abbastanza basso (per problemi legati al costo ed allo spazio). di C.Bettini, A.Guiggiani e T.Lorenzetti 113 11.2 Convertitori Figura 11.8.: Modulatore sigma-delta (Σ∆) Figura 11.9.: Struttura di un filtro FIR Convertitori subranging sono convertitori (Figura 11.11) che hanno in ingresso un Sample & Hold (S&H) che diminuisce notevolmente la banda passante (e quindi la velocità di conversione). Di base i convertitori subranging sono composti da due convertitori flash tali che il primo codifica in bassa risoluzione ed il secondo effettua un secondo passaggio in modo tale da avere risoluzione maggiore. In questo modo con due convertitori flash da 4 bit (24 = 16 comparatori) si riesce a codificare un numero a 8 bit (che altrimenti necessiterebbe di 28 = 256 comparatori). Convertitori pipeline sono (Figura 11.12) formati da una cascata di convertitori subranging separati da Sample & Hold. Prestazioni di un ADC le prestazioni di un convertitore analogico-digitale sono condizionate principalmente da: 1. rumore dei circuiti analogici, filtrabile con un filtro passa-basso; 114 di C.Bettini, A.Guiggiani e T.Lorenzetti 11.2 Convertitori Figura 11.10.: Convertitore Flash 2. rumore dell’ADC, non filtrabile ma comunque limitabile con accorgimenti di progetto (l’equivalente del “noise shaping” per il converitore Σ∆); 3. jitter di apertura; 4. non linearità dei circuiti analogici; 5. qualità del clock. Dato che per un generico convertitore a N bit si ha un rapporto segnale/rumore: SN R = 6N + 1.8 è possibile definire il numero effettivo di bit (ENOB, quello che poi si può usare come riferimento per le specifiche) come: EN OB = SIN AD − 1.8 6 di C.Bettini, A.Guiggiani e T.Lorenzetti 115 11.2 Convertitori Figura 11.11.: Convertitore Subranging dove il SIN AD (Signal-to-Noise And Distortion Ratio) è il parametro che rappresenta la qualità del segnale uscente da un dispositivo. Solitamente quando si deve scegliere il numero di bit che il nostro convertitore deve avere per soddisfare le specifiche si considera un +2 sulle specifiche stesse (quindi se sono richiesti 10 bit si usa un convertitore a 12 bit). 116 di C.Bettini, A.Guiggiani e T.Lorenzetti 11.2 Convertitori Figura 11.12.: Convertitore Pipeline di C.Bettini, A.Guiggiani e T.Lorenzetti 117 12. Conversione Digitale-Analogico (DAC) Prima di studiare la conversione digitale-analogico è necessario introdurre il funzionamento e le caratteristiche in frequenza degli interruttori. Interruttori analogici NMOS i primi interruttori che analizzeremo sono quelli basati su un NMOS (Figura 12.1). Questi interruttori sono caratterizzati da una relazione ingresso-uscita: Vout = Rs Rs + Rl + Req dove Req è la resistenza equivalente del transistor NMOS caratterizzata da un andamento come quello in Figura 12.2: parte da un valore in condizione di ON i−1 h Ron = 2k W (V − V ) dove k, W, L sono costanti e dimensioni caratteristiche dd t L del transistor, e tende (in condizione di OFF) ad una Rof f = 1 2k W L (Vgs − Vt ) Ricordiamoci di prestare attenzione al fatto che nell’istante di switching Req (e, di conseguenza, Vout ) non ha alcun comportamento lineare. Interruttori analogici CMOS analizziamo adesso gli interruttori CMOS (Figura 12.3): questi non sono altro che un NMOS accoppiato ad un PMOS e, di conseguenza, la resistenza equivalente è data (Figura 12.4) da quella di un NMOS unita a quella 119 Conversione Digitale-Analogico (DAC) Figura 12.1.: Interruttore basato su NMOS Figura 12.2.: Resistenza equivalente per interruttori NMOS di un CMOS (che è un NMOS “ribaltato”). Si ottiene così un andamento come in Figura 12.4 in cui al variare di V la Req è più o meno costante. Inoltre un interruttore CMOS essendo un parallelo tra un NMOS ed un PMOS è caratterizzato da una resistenza equivalente che è il parallelo tra due resistenze simili (e quindi all’incirca la metà). Analisi in frequenza andiamo ora ad analizzare in frequenza il comportamento di un interruttore. Iniziamo dalla condizione di ON, Figura 12.5.a Notiamo subito che saranno presenti due capacità di accoppiamento: una Cin tra source e terra, l’altra Cout tra terrae drain. In generale Cout è maggiore e, quindi, 120 di C.Bettini, A.Guiggiani e T.Lorenzetti Conversione Digitale-Analogico (DAC) Figura 12.3.: CMOS Figura 12.4.: Resistenza equivalente per interruttori CMOS comanda la frequenza di taglio in condizione di ON (Figura 12.6.a) è pari a: ωON = 1 Req Cout Al contrario in condizione di OFF (Figura 12.5.b) prevale il comportamento capacitivio dell’interruttore, quindi si ottiene una Cof f molto piccola che porta ad avere una frequenza di taglio in condizione di OFF (Figura 12.6.b) definita come: ωOF F ≈ 1 ωON 10 Si determina così la risposta in frequenza di un interruttore che si comporta da filtro di C.Bettini, A.Guiggiani e T.Lorenzetti 121 Conversione Digitale-Analogico (DAC) passa-basso in fase di ON e in filtro passa-alto in fase di OFF. Figura 12.5.: Circuito equivalente per interruttori Figura 12.6.: Risposta in frequenza di un interruttore Test funzionamento dinamico ADC Per il test dinamico di un ADC si può misurare il rapporto tra Potenza del Segnale, Ps , e Potenza di Rumore, PN , valutando separatamente i 2 contributi nello spettro. Per fare ciò (in frequenza) occorre fare riferimento ad un segnale sinusoidale. Tale segnale deve essere: • coerente, quindi fin = fsegnale k . Il risultato può però dipendere da k; • far rientrare un numero primo di cicli nella finestra di FFT. Se non si dispone di sorgenti con le necessarie accuratezza e stabilità in frequenza si può: • per evitare i problemi del campionamento coerente, si possono usare frequenze “incoerenti” aggiungendo una finestra di pesaggio che abbassi i lobi laterali; • individuare quanti punti contribuiscono alla potenza del segnale, in base all’SNR previsto. 122 di C.Bettini, A.Guiggiani e T.Lorenzetti 12.1 Convertitori DAC e sintesi di segnali 12.1. Convertitori DAC e sintesi di segnali Data una Vref di riferimento, un convertitore DAC segue lo schema in Figura 12.7. Principalmente questa tecnica è soggetta a quattro limitazioni: 1. la tensione Vref ha un limite sulla massima frequenza poichè gli switch hanno una caratteristica capacitiva che li fa lavorare in frequenza come dei filtri passa-basso; 2. lo slew-rate massimo è legato al settling-time, che dipende a sua volta sia dal DAC che dalla rete analogica (AO) a valle del convertitore. Solitamente si definisce il settling-time ts come: ts = q t2sDAC + t2sAO 3. è possibile che si abbiano degli spike sull’uscita (detti “glitch” di corrente) in corrispondenza di variazioni massive di bit (ad esempio da 0001 a 1110); questi sono dovuti principalmente al fatto che gli interruttori al momento della commutazione (simultanea soprattutto) immettono dei piccoli spike di corrente sul circuito che si sommano e influenzano il comportamento del convertitore; 4. è presente inoltre un fenomeno di frequency-roll-off dovuto alla natura di tipo S&H del DAC. Questo porta ad avere variazioni di ampiezza del segnale in funzione della frequenza, ed è possibile risolverlo con tecniche di pre-enfasi e di post-moltiplicazione. Sintesi di segnali andiamo adesso ad analizzare come è possibile sintetizzare un segnale usando un DAC. Se utilizziamo una struttura come quella in Figura 12.8 è semplice ricavare la frequenza massima del segnale in uscita come: fout = fck 2N quindi per variare fout è necessario variare fck . Per sintetizzare un segnale digitale è possibile utilizzare una configurazione come quella in Figura 12.9. Tale configurazione è caratteristica del DDS (Direct Digital di C.Bettini, A.Guiggiani e T.Lorenzetti 123 12.1 Convertitori DAC e sintesi di segnali Figura 12.7.: Principio di funzionamento di un DAC Synthesizer) che è caratterizzato dal fatto che: fout = K fck 2N quindi per variare fout basta variare l’indice K. Per sintetizzare un segnale con un DSS si procede come segue: 1. si ipotizza il valore massimo di fck compatibile con la temporizzazione; 2. si valuta 2N = fck ; ∆f se risulta eccessivo, ridurre fck ; 3. si valuta il rapporto tra frequenze armoniche e fondamentali; 4. si valuta il rapporto tra ampiezze fondamentali e armoniche ottenibile col filtro di ordine assegnato. Se eccessivo, aumentare Fck. Facciamo un esempio: data una EPROM da 8Kx8, ovvero con N = 13 (213 = 8192 = 8K) ed m = 8 (deriva dal x8 dato che la EPROM ha in ingresso n = 13 ed in uscita m = 8), si desidera fout ∈ [1 ÷ 100] KHz, tAA = 50ns ed un SNR>50dB. Dalla definizione della fout si ricava che ∆f 6 1KHz da cui: ∆f = 124 fck 6 1KHz ⇒ fck 6 8M Hz = fckM AX 2m di C.Bettini, A.Guiggiani e T.Lorenzetti 12.1 Convertitori DAC e sintesi di segnali Figura 12.8.: Sintesi di segnali tramite DAC Figura 12.9.: Sintesi di segnali avanzata tramite DAC Ipotizziamo ora di applicare un filtro LPF del terzo ordine (che guadagna quindi 18dB per raddoppio di frequenza) con una frequenza di taglio a 500KHz. Avremo: 500KHz → 1M Hz18dB → 2M Hz36dB → 4M Hz54dB quindi per fck = 4M Hz (che è una frequenza minore di fckM AX ) si soddisfa il requisito sull’SNR. Andiamo adesso ad analizzare il vincolo su: fckM AX 6 1 tckmin dove tckmin è dato da una serie di contributi: tckmin = tcq + tsu + tAA di C.Bettini, A.Guiggiani e T.Lorenzetti 125 12.1 Convertitori DAC e sintesi di segnali dove tcq e tsu è plausibile ipotizzarli pari a 10ns e tAA è una specifica nota. Quindi tckmin = 70ns da cui deriva il vincolo: fckM AX 6 12.5M Hz Modulazione e sintesi di segnali in quest’ultimo esempio (Figura 12.10) si vanno ad analizzare le tecniche di modulazione per segnali sintetizzati tramite DAC. È infatti possibile applicare una modulazione in frequenza (FM), una modulazione in fase (PM) oppure una in ampiezza (AM) a seconda di dove si fanno entrare determinati ingressi (additivi, nel caso FM e PM, e moltiplicativi, nel caso AM). Figura 12.10.: Sintesi di segnali modulata tramite DAC 126 di C.Bettini, A.Guiggiani e T.Lorenzetti 13. Metastabilità e distribuzione di clock La metastabilità riguarda tutti quei problemi a livello di registro dovuti ad asincronie nel sistema digitale. Di norma, infatti, il registro funziona se vengono rispettati i tempi di setup e di hold durante i quali i dati devono essere disponibili. Come si può notare in Figura 13.1, nell’istante in cui si accede ad un dato nel registro è necessario tenere conto di un ritardo td (a sua volta comprensivo del tempo di setup tsu , l’intervallo piccolo in figura) e di un tempo di hold th . Più precisamente, la somma tsu + th determina una zona proibita entro la quale il dato non può (o, meglio, non deve) interrompersi. La metastabilità si inizia a notare se facciamo diminuire td fino ad entrare nella zona proibita: è così possibile determinare un tempo critico tale che, per td < tcrit , la variazione dell’ingresso viene ignorata in uscita con conseguente piena metastabilità. Si noti che solitamente la zona proibita è nell’ordine del nanosecondo ed è specificata nei datasheet dei componenti; è possibile inoltre ricavare (sperimentalmente) la zona critica (entro la quale siamo certi di avere metastabilità): questa è solitamente nell’ordine dei 100 ps e non viene fornita dal costruttore all’interno del datasheet. Figura 13.1.: Caratteristiche temporali della metastabilità 127 Metastabilità e distribuzione di clock Infatti la metastabilità può manifestarsi in diversi modi (Figura 13.2). Più precisamente si nota: • l’uscita non affetta da metastabilità (Figura 13.2.1); • una leggera perdita del bordo quando ci si avvicina allo stato metastabile (Figura 13.2.2); • una perdita più consistente del bordo avvicinandoci ancora alla criticità (Figura 13.2.3); • un primo stadio in cui si hanno delle failure (rare) alternate a frequenti stati ben rilevati (Figura 13.2.4); • lo stato metastabile caratterizzato da failure frequenti (Figura 13.2.5); • metastabilità piena: nessuno stato viene rilevato (Figura 13.2.6). Figura 13.2.: Fasi della metastabilità MTBF e metastabilità per ogni componente è possibile definire un indice di MTBF (Maximum Time Between Failure) dato da: M T BF = eT ·td fck · fin · T0 dove T, T0 sono caratteristici dell’interruttore. Si riporta in Figura 13.3 l’andamento dell’MTBF per vari componenti al variare del parametro td (su scala logaritmica). Normalmente le logiche più rapide e ad alta 128 di C.Bettini, A.Guiggiani e T.Lorenzetti Metastabilità e distribuzione di clock tecnologia sono quelle con pendenza più elevata (nel grafico sulla sinistra) mentre quelle standard sono quelle più “piane” (si noti la 74Standard). Figura 13.3.: Esempio di MTBF per vari componenti al variare di td Rimedi alla metastabilità un possibile rimedio alla metastabilità è quello di utilizzare un circuito a doppia sincronizzazione (Figura 13.4). In questo modo non si ha metastabilità in quanto la risposta di QF F 2 è decisa e non metastabile (però può essere errata, e qua è un po’ difficile metterci una pezza): in parole povere non si hanno stati metastabili ma si può comunque incorrere in errori. Definite le frequenze fin1 ed fin2 come in figura, il valore dell’MTBF diventa: M T BF (2) = eT ·td fin2 · fck · T0 dove: fin2 = 1 fin1 · T0 · fck = T 1 M T BF (1) e fck di C.Bettini, A.Guiggiani e T.Lorenzetti 129 Metastabilità e distribuzione di clock quindi la relazione tra M T BF (2) (del doppio sincronizzatore) ed M T BF (1) (della logica singola) è: M T BF (2) ≈ 3M T BF (1) Figura 13.4.: Circuito a doppia sincronizzazione Rilevare la metastabilità per rilevare la metastabilità è possibile porre il DUT (Device Under Test) in un circuito come quello in Figura 13.5: se l’uscita è 0 non vi è metastabilità; quando l’uscita diventa 1 viene rilevato uno stato metastabile. In questo modo è possibile contare gli stati metastabili e ricavare l’MTBF della logica. Figura 13.5.: Circuito a utilizzato per rilevare la metastabilità 130 di C.Bettini, A.Guiggiani e T.Lorenzetti 13.1 Distribuzioni di clock 13.1. Distribuzioni di clock Tipicamente il clock deve essere distribuito a diversi ricevitori dello stesso sistema ed arrivare a tutti simultaneamente. Si definisce così “clock skew” l’incertezza nel tempo di arrivo del clock. Possibili cause dello skew sono le tolleranze delle logiche, la temperatura di funzionamento e la tensione di alimentazione dal punto di vista fisico. Per quanto riguarda le cause legate alle interconnessioni si elencano la differente lunghezza dei percorsi di clock ed il diverso ritardo di propagazione nonchè le variazioni del tipo di linea (strisce, microstrisce,...) o delle costanti dielettriche. Alcune tecniche per sopperire a questi problemi sono: 1. utilizzo di buffer: questa tecnica però porta con sè le problematiche legate ai buffer; i buffer, infatti, anche se uguali (almeno su carta) sono caratterizzati da ritardi output-to-output differenti che introducono comunque uno skew; 2. allungamento delle piste, cosa fattibile ma per piccoli ritardi (Figura 13.6); 3. disposizione di clock “a stella” (Figura 13.7) rispetto ad un generatore di clock centrale al circuito; 4. utilizzo di generatori di clock (Figura 13.8) in aggiunta ai quali può comunque essere utile utilizzare la tecnica di allungamento delle piste. Figura 13.6.: Esempio di allungamento variabile delle piste Generazione di un clock stabile il problema fondamentale resta comunque quello di riuscire a generare un clock stabile e la cui frequenza sia variabile a piacere. Per di C.Bettini, A.Guiggiani e T.Lorenzetti 131 13.1 Distribuzioni di clock Figura 13.7.: Esempio di allungamento delle piste e distribuzione di clock a stella fare ciò si ricorre al PLL (Phase Locked Loop), Figura 13.9. Si noti che a valle del PLL è presente un “programmable delay” che varia lo skew dei segnali di clock generati. Analizzando un PLL internamente (Figura 13.10) si nota che la frequenza di uscita Fout è data da: Fout = P · Fin Q Tutto ciò è realizzato grazie ad una retroazione che assicura la stabilità del clock e ad un filtro passa basso (LPF) che assicura un clock pulito. Si noti in Figura 13.12 una possibile realizzazione fisica del phase-detector che altro non fa se non identificare lo sfasamento tra i due segnali A e B. Infine, in Figura 13.12, si mostra una possibile realizzazione fisica del VCO, l’oscillatore controllato dalla tensione Vc e che va a generare l’uscita out. Si riporta in Figura 13.13 un esempio di “programmable delays”: mentre il PLL genera un clock stabile, questo secondo componente deve applicare uno skew alle varie fck del generatore. Il modello Roboclock (in Figura 13.13) fu sviluppato alla fine degli anni ’80 ed è un dispositivo multi-output che attraverso dei buffer ritarda i 132 di C.Bettini, A.Guiggiani e T.Lorenzetti 13.1 Distribuzioni di clock Figura 13.8.: Esempio di clock-generator Figura 13.9.: Esempio di PLL segnali di clock (generando un clock-skew programmabile che ci permette di risolvere i problemi di sincronizzazione). di C.Bettini, A.Guiggiani e T.Lorenzetti 133 13.1 Distribuzioni di clock Figura 13.10.: Interno di un PLL Figura 13.11.: Realizzazione di un phase-detector Figura 13.12.: Realizzazione del VCO 134 di C.Bettini, A.Guiggiani e T.Lorenzetti 13.1 Distribuzioni di clock Figura 13.13.: Roboclock di C.Bettini, A.Guiggiani e T.Lorenzetti 135 14. Sistemi ad alta velocità Andiamo adesso ad introdurre i sistemi ad alta velocità. Analizziamo prima di tutto la risposta in frequenza di un segnale digitale, e notiamo che questo si comporta come un filtro passa-basso (Figura 14.1) in cui: f1 = 1 1 e f2 = πtw πtr Figura 14.1.: Analisi in frequenza di segnali digitali Linee di trasmissione analizziamo adesso le linee di trasmissione più utilizzate e le loro caratteristiche: 137 Sistemi ad alta velocità 1. microstriscia (Figura 14.2.a): è tale per cui valgono le relazioni 87 5.96h Z0 = √ ln 0.8w + t r + 1.41 ! ≈ 50Ω e: √ τp = 1.016 0.475r + 0.67 [ns/f t] 2. strip-line (Figura 14.2.b): in questo caso valgono 60 4b Z0 = ln r 0.67π (0.8w + t) ! 3. bus-line (Figura 14.3): sono da analizzare in modo un po’ più accurato. Figura 14.2.: Microstriscia (a) e strip-line (b) Figura 14.3.: Esempio di Bus-line Analizziamo adesso cosa succede su una generica linea di trasmissione (Figura 14.4) dal punto A al punto B, considerata la resistenza interna del generatore (solitamente 138 di C.Bettini, A.Guiggiani e T.Lorenzetti Sistemi ad alta velocità RG ≈ 50Ω) e quella del carico (RL ) dipendente dal tipo di quest’ultimo: se il segnale entra in un MOS, ad esempio, la sua resistenza interna è assimilabile ad RL = ∞. Nel tragitto A → B l’onda parte da A in t = 0 ed arriva in B nell’istante t = τp . A questo punto da B parte un’inda riflessa verso A caratterizzata da un coefficiente di riflessione: ρB = RL − Z0 RL + Z0 dove Z0 è l’impedenza equivalente della linea. Allo stesso modo nell’istante t = 2τp l’onda riflessa di B arriva in A producendo una seconda onda riflessa (verso B) caratterizzata da un coefficiente di riflessione: ρA = RG − Z0 RG + Z0 Figura 14.4.: Generalizzazione di una linea di trasmissione Quindi in ogni istante rispettivamente per t = 0, τp , 2τp , 3τp , ... su ogni estremo (rispettivamente A e B) sono presenti (Figura 14.6): • onda incidente; • onda riflessa; • onda pre-esistente; il tutto partendo da determinate condizioni iniziali in cui si considera che la tensione VA è pari a: VA |t=0 = VG Z0 RG + Z0 di C.Bettini, A.Guiggiani e T.Lorenzetti 139 Sistemi ad alta velocità Figura 14.5.: Esempio di analisi di linee di trasmissione e dalle condizioni a regime in cui, non essendoci più riflessioni, scompare il termine Z0 e si ottiene: VA |t=∞ = VB |t=∞ = V̄ = VG RL RG + RL Si può così utilizzare il metodo grafico di Bergeron calcolando i punti di partenza VA (0) e di arrivo VA (∞) e facendo “rimbalzare” l’onda (con coefficiente angolare Z0 ) sulla curva del generatore (caratterizzata da un coefficiente angolare −RG ) e sulla curva del carico (con coefficiente angolare RL ) fino a raggiungere V̄ (un esempio in Figura 14.6). Si noti che se il carico è una logica (RL = ∞) la retta del carico coincide con l’asse y mentre se il carico è un cortocircuito la retta di carico coincide con l’asse x (ed ha poco senso). Naturalmente a seconda del rapporto in cui sono RG , Z0 ed RL avremo diversi grafici. Analizziamo in Figura 14.7 il caso in cui RL = ∞ (quindi plausibilmente sul 140 di C.Bettini, A.Guiggiani e T.Lorenzetti 14.1 Terminazione di linee digitali Figura 14.6.: Esempio di analisi con i diagrammi di Bergeron carico abbiamo delle porte logiche) e in Figura 14.8 altri casi interessanti. Andiamo adesso ad analizzare un metodo per la misura dell’impedenza di linea (Figura 14.9): data R0 = 25 in quanto pari al parallelo tra la resistenza dell’oscilloscopio (solitamente 50Ω) e del generatore di impulsi (50Ω), si può calcolare: Z0 = R0 −1 V0 VA estrapolando i dati mancanti dal grafico, e τ · Z0 = L e τ Z0 = C. 14.1. Terminazione di linee digitali Le riflessioni si possono combattere terminando una linea, azione che può essere compiuta in parallelo (al termine della linea) od in serie (all’inizio della linea). Le terminazioni di tipo parallelo sono tendenzialmente più “power consuming”; proprio per questo le terminazioni di tipo serie sono preferite nella maggior parte dei casi, eccetto che nel controllo di bus (poichè con queste terminazioni il circuito è altamente corrotto da fenomeni di metastabilità). Consideriamo la generica porta logica (Figura 14.11) e analizziamo gli effetti che una terminazione avrebbe sull’uscita: di C.Bettini, A.Guiggiani e T.Lorenzetti 141 14.1 Terminazione di linee digitali Figura 14.7.: Dipendenza tra RG e Z0 con RL = ∞ • una terminazione serie (adattando con R = Z0 + RG , Figura 14.10) è generalmente preferita nelle linee porta-porta ma non può essere applicata ai bus poichè porta alla metastabilità; • una terminazione parallelo come quella in Figura 14.12 ha come contro un alto consumo di potenza; inoltre, considerata la struttura della logica (Figura 14.11), una terminazione di questo tipo fa si che i due transistor non lavorino allo stesso modo (in pratica il transistor alto lavora di più e, quindi, si scalda). Per sopperire a questo problema è possibile utilizzare una terminazione con R = 2Z0 , così da avere un mismatch del 50% (accettabile) ed un minor dispendio di potenza (comunque i transistor della logica continuano a lavorare in modo errato); • una terminazione parallelo con partitore (detta rete di Thevenin, Figura 14.13), se è tale che il parallelo tra R1 ed R2 è circa 2Z0 (con R2 > R1 così da essere meno dispersivo) allora la logica lavora in modo simmetrico. Resta comunque il problema delle perdite di potenza; • una terminazione parallelo con carico R-C (Figura 14.14) è tale per cui non si rischia più la metastabilità poichè a regime il condensatore equivale ad un gene- 142 di C.Bettini, A.Guiggiani e T.Lorenzetti 14.1 Terminazione di linee digitali Figura 14.8.: Dipendenza tra RG , RL e Z0 ratore di tensione (di V2 ). Sulla resistenza scorrono gli altri V2 così da dimezzare la potenza dissipata (anche per questo è uno dei metodi di terminazione più usati). Vige la regola pratica secondo cui R · C = 4τp . n di C.Bettini, A.Guiggiani e T.Lorenzetti 143 14.1 Terminazione di linee digitali Figura 14.9.: Esempio di misura dell’impedenza di linea Figura 14.10.: Terminazione serie Figura 14.11.: Porta logica 144 di C.Bettini, A.Guiggiani e T.Lorenzetti 14.1 Terminazione di linee digitali Figura 14.12.: Terminazione parallelo 1 Figura 14.13.: Terminazione parallelo con partitore Figura 14.14.: Terminazione parallelo R-C di C.Bettini, A.Guiggiani e T.Lorenzetti 145 15. Rumori e circuiti ad alta velocità Andiamo adesso ad analizzare le problematiche relative ai circuiti ad alta velocità. Nello specifico, cerchiamo di studiare i disturbi di commutazione (switching noise), l’accoppiamento tra le linee (crosstalk) e gli approcci al “buon progetto” di una scheda a circuiti stampati (PCB). 15.1. Switching-noise Lo switching noise è un fenomeno di scarica che avviene nella commutazione dei transistor interni alle logiche nel passaggio alto-basso. Se ciò avviene in modo non perfettamente sincrono, infatti, è possibile che si generi un glitch (picco) in corrispondenza della commutazione (Figura 15.1). Se consideriamo il fatto che l’uscita, all’istante di commutazione, si scarica su una capacità parassita C tale che: ∆i = C dV dt ad esempio per C = 25pF , V = 5V e dt = trise = 2ns si ottiene un ∆i ≈ 60mA che comunque va considerato come un valore medio. Essendo un valore medio, se consideriamo che uno spike di corrente ha una forma d’onda triangolare, ne deriva che il valore massimo di tale spike (imax ) è circa 120mA (non poco). Se adesso si considerano le induttanze parassite della linea, si ottiene una relazione del tipo: ∆V = L di dt 147 15.1 Switching-noise che aumenta all’aumentare della velocità della logica (quindi della frequenza di switching) per diminuire con il “power consumption” (la potenza della logica). Se riprendiamo i 120mA calcolati in precedenza, considerata una induttanza media della linea (L0 in Figura 15.2) di 5nH/cm si ottiene un ∆V compreso tra 0.5V ed 1V che, su una logica che lavora a 5V , è un disturbo assolutamente non trascurabile (che fra l’altro va a spostare la tensione di riferimento a massa, Figura 15.3, e può portare a fenomeni di metastabilità sulle “uscite quiescenti”). Alcune soluzioni per prevenire lo switching noise sono: 1. aumentare la velocità e, quindi, la tecnologia della logica (ma non troppo poichè come visto aumentando la velocità cresce ∆V ); 2. ottimizzare il progetto in modo da avere bus non flottanti e pin non disconnessi; 3. non utilizzare i “socket” che vanno ad aumentare la capacità parassita; 4. progettare in modo corretto i piani di alimentazione e di massa in modo da minimizzare gli accoppiamenti; 5. limitare quanto più possibile la capacità del carico. Figura 15.1.: Analisi dello switching noise 148 di C.Bettini, A.Guiggiani e T.Lorenzetti 15.2 Crosstalk Figura 15.2.: Rappresentazione della linea con la sua induttanza caratteristica (L0 ) Figura 15.3.: Esempio di switching noise accentuato dalle commutazioni simulanee di più logiche 15.2. Crosstalk un altro fenomeno che, assieme allo switching noise, caratterizza i problemi dei circuiti ad alta velocità, è l’accoppiamento tra le piste (o crosstalk). Questo fenomeno, come il precedente, è prevenibile in fase di progetto adottando buone tecniche di progetto del layout delle schede digitali. Il crosstalk avviene quando una linea di trasmissione (detta “linea killer”) si accoppia ad una seconda (“linea vittima”) con la quale non dovrebbe avere niente a che fare. Ciò avviene in quanto una linea di trasmissione è sempre rappresentabile a parametri concentrati attraverso un insieme di induttanze e capacità di accoppiamento (Figura 15.4 e Figura 15.7). Il crosstalk può essere di due tipi a seconda di chi è l’ingresso (killer) e di chi è di C.Bettini, A.Guiggiani e T.Lorenzetti 149 15.2 Crosstalk l’uscita (vittima): • forward-crosstalk (Figura 15.5): in questo caso vi è un ritardo fisso (pari a τp ) in cui arrivano i disturbi indipendentemente dal percorso. Ne derivano così due componenti (una capacitiva ed una induttiva) che si sommano annullandosi tra loro (data la loro natura, Figura 15.7) se i componenti sono omogenei; • reverse-crosstalk (Figura 15.6): in questo caso il problema è più complesso in quanto i disturbi arrivano in istanti che dipendono dal cammino intrapreso, quindi vanno a sommarsi tra loro introducendo un bel disturbo. Ma come si misura il crosstalk? Idealmente è misurabile facendo il rapporto tra segnale di disturbo sulla linea vittima e segnale sulla linea killer. In realtà ci viene in aiuto una formula empirica: vnoise k = 2 D ∆vkiller 1− H dove D è la distanza tra killer e vittima, H la distanza tra vittima e massa e k caratteristica della struttura della linea. Se andiamo ad analizzare il reverse-crosstalk (Figura 15.8) si va a definire una lunghezza critica della linea Lc tale che: • per L < Lc (Figura 15.8.a): solo pochi impulsi (di durata tr ) contribuiscono all’integrale (di durata tr + 2τp ) (tutti gli impulsi che contribuiscono sono parzialmente sovrapposti); • per L = Lc (Figura 15.8.b): l’integrale raggiunge il valore massimo e subito dopo inizia a decrescere (tutti gli impulsi che contribuiscono sono parzialmente sovrapposti, eccetto l’ultimo); • per L > Lc (Figura 15.8.c): l’integrale raggiunge il valore massimo e lo mantiene per 2τp − tr (non tutti gli impulsi che contribuiscono sono parzialmente sovrapposti). Si noti che il reverse crosstalk è massimo in questo caso e, al contempo, che il caso L > Lc impone anche che tr < 2τp (cosa da tenere in conto al momento della progettazione). Per quanto riguarda i metodi di prevenzione per il crosstalk si osserva che: 150 di C.Bettini, A.Guiggiani e T.Lorenzetti 15.2 Crosstalk 1. il piano di massa deve essere il più integro possibile; 2. è bene utilizzare delle linee di guardia per minimizzare i fenomeni di accoppiamento; 3. è bene minimizzare gli accoppiamenti tra la linea che porta il clock (il killer per antonomasia) e le altre linee; 4. è bene non porre ingressi ed uscite vicini (così da limitare i casi di reverse crosstalk). Figura 15.4.: Rappresentazione di una linea di trasmissione a parametri concentrati Figura 15.5.: Forward crosstalk di C.Bettini, A.Guiggiani e T.Lorenzetti 151 15.3 Problematiche di Layout Figura 15.6.: Reverse crosstalk Figura 15.7.: Caratteristiche del crosstalk 15.3. Problematiche di Layout Analizziamo adesso il fatto che, ad ogni commutazione, si possono avere elevati di valori di dt e di dv i quali„ coinvolgendo induttanze e capacità parassite, possono dt di dar luogo a notevoli fenomeni di disturbo. Un modo approssimato per stimare dt , dove dt = tr , è il seguente. Considerato che: ∆i = C 152 dv dt di C.Bettini, A.Guiggiani e T.Lorenzetti 15.3 Problematiche di Layout Figura 15.8.: Analisi del reverse crosstalk possiamo sostituire in: ∆v = L dove d2 v dt2 di ∆i d2 v =L =L·C · 2 dt tr dt ha un andamento come in Figura 15.9. Si noti che è possibile approssimare: ( di max dt ) d2 v = C · max dt2 ( ) = 1.5 · C · ∆V t2r che per dei valori plausibili di C = 50pF , ∆V = 3.7V e tr = 2ns porta a picchi massimi di 70mA/ns che non possiamo ignorare. Figura 15.9.: Analisi dell’andamento di un disturbo ∆v Ci si chiede adesso quale sia il percorso di “richiusura” delle correnti di disturbi. Si individuano due casi: 1. se la natura del circuito è di tipo DC (o comunque a frequenze molto basse) di C.Bettini, A.Guiggiani e T.Lorenzetti 153 15.3 Problematiche di Layout predomina la natura resistiva della linea ed il loop si richiude sul percorso a resistenza minima (quindi sul più breve); 2. se la natural del circuito è di tipo AC (e ciò avviene sempre nel caso di circuiti ad alta frequenza) predominano le caratteristiche induttive sulla linea ed il loop si richiude sul percorso a minima induttanza (e, quindi, ad area minima). Proprio per questo è dannoso togliere integrità al piano di massa: in tal modo vi è la possibilità che le correnti si richiudano al di sotto del piano di massa. Notiamo ora che l’induttanza parassita di un circuito aumenta con la lunghezza e diminuisce con la sezione, limitandosi nel caso in cui sia presente anche un piano di massa. Per porre un freno a queste problematiche è possibile applicare una capacità di bypass (Cb , in parallelo alle logiche) in modo tale da ridurre l’induttanza parassita andando a limitare l’area del loop (un esempio applicativo è quello di utilizzare una PCB a quattro strati due dei quali sono massa ed alimentazione; in questo modo si ha all’incirca 100pF/inch). Riprendendo: i(t) = C · dv dt e considerando dt = tr ed i(t) = i0 costante, si ottiene la condizione: Cb > Cmin = i0 · tr dvmax Allo stesso modo è possibile applicare, in serie alla linea (a monte quindi), una induttanza di blocco che va a filtrare il rumore uscente (verso l’alimentazione) e parte di quello entrante (alle alte frequenze). Di conseguenza sono possibili soluzioni di progetto: 1. mantenere il piano di massa il più integro possibile; 2. controllare in fase di progetto i percorsi di ritorno; 3. seguire le regole per ridurre switching noise e crosstalk (che comunque influenzano in modo pesante il funzionamento della scheda); 4. usare più Cb piccoli anzichè un unico Cb di elevato valore; 154 di C.Bettini, A.Guiggiani e T.Lorenzetti 15.3 Problematiche di Layout 5. separare le alimentazioni applicando delle induttanze di blocco. di C.Bettini, A.Guiggiani e T.Lorenzetti 155 Parte IV. Relazione di macchine elettriche 157 15.4 Introduzione e Misure 15.4. Introduzione e Misure 15.5. Obiettivi L’obiettivo della presente relazione è la caratterizzazione di un motore asincrono trifase e la costruzione del relativo diagramma circolare. Per poter ricavare i dati necessari, bisogna eseguire tre prove sulla macchina: • Misura delle resistenze statoriche; • Prova a vuoto (rotore senza carico); • Prova in corto circuito (rotore bloccato). 15.6. Dati di targa del motore trifase Pn Vn 5,5 kW 380 V In cos ϕ η n1 classe d’isolamento 11,7 A 0,83 88,3 % 1500 B Tabella 15.2.: Dati di targa 15.7. Misura della resistenza statorica Per determinare la resistenza di fase statorica si esegue una prova in corrente continua andando ad usare l’inserzione volt-amperometrica con il voltmetro a valle. La corrente applicata deve essere inferiore al 10% della corrente nominale al fine di ridurre al minimo il riscaldamento degli avvolgimenti. Si determina così la resistenza delle fasi concatenate come media delle tre resistenze rilevate sulle tre coppie di morsetti. La resistenza di fase, essendo i morsetti collegati a stella, viene infine ottenuta come Rmedia/2. di C.Bettini, A.Guiggiani e T.Lorenzetti 159 15.8 Tabelle delle misure 15.8. Tabelle delle misure Misure delle resistenze morsetti I [A] V [V] R [Ω] ∆R% 1-2 0,98 1,27 1,30 ±0, 3% 1-3 0,99 1,29 1,30 ±0, 3% 2-3 0,98 1,27 1,30 ±0, 3% Tabella 15.4.: Misure effettuate 15.9. Dati ricavati dalle misure Nel nostro caso i valori ottenuti sono stati: Rmedia = Rstat = 1, 30 + 1, 30 + 1, 30 = 1, 30Ω ± 0, 3% 3 Rmedia = 0, 65Ω ± 0, 3% 2 Si riporata quindi il valore della Rstat alla temperatura di riferimento convenzionale Tr = 75°dipendente dalla classe d’isolamento della macchina (B). Il fattore di riporto in temperatura si ottiene dalla seguente formula: KT = 234, 5 + Tr 234, 5 + Tp quindi la resistenza riportata alla temperatura convenzionale è: T Rstat = K T Rstat 160 di C.Bettini, A.Guiggiani e T.Lorenzetti (15.1) 15.9 Dati ricavati dalle misure I risultati sono riportati in Tabella 15.6. Temperatura di prova Tp 21°C ±0, 5% Temperatura convenzionale Tr 75°C Fattore di temperatura KT 1,21 ±0, 5% T Resistenza statorica a Tr Rstat 0,79 ±0, 8% Tabella 15.6.: Resistenza riportata alla temperatura convenzionale di C.Bettini, A.Guiggiani e T.Lorenzetti 161 16. Prova a vuoto La prova a vuoto consiste nel far funzionare il motore alla tensione nominale senza applicare nessun carico al rotore, quindi misurare la potenza a vuoto (P0 ), la corrente a vuoto (I0 ) e il fattore di potenza (cos ϕ0 ). La prova viene effettuata inizialmente alimentando il motore ad un voltaggio pari al 110% della tensione nominale; a questo punto si va a decrementare la tensione applicata al motore fino ad arrivare al 25% circa di Vn . 16.1. Tabelle delle misure prova n° V0 [V] I0 [A] P0 [W] cos ϕ0 1 2 398,6 381 5,328 4,796 451,3 403,2 0,123 0,127 Tabella 3 4 5 6 7 8 9 10 343,7 297,1 264,7 221,3 182,5 159,7 138 121,2 3,961 3,21 2,759 2,241 1,855 1,628 1,486 1,357 334,8 275,6 243,1 202,5 175,6 158,1 142,3 130,3 0,142 0,167 0,192 0,236 0,299 0,351 0,401 0,457 16.2.: Dati ricavati dalle misure a vuoto 16.2. Dati ricavati dalla prova a vuoto Nella prova a vuoto le perdite nel ferro e nel rame rotoriche possono essere trascurate poichè le correnti indotte nel rotore sono molto ridotte. Quindi: P0 ' Pj0 + Pf m = 403.2 ± 1, 1% 163 16.2 Dati ricavati dalla prova a vuoto Ovvero la potenza assorbita a vuoto è la somma delle potenze perse nel rame (e nel ferro) statorico e delle potenze meccaniche. Le perdite nel rame statorico possono essere espresse come: 3 Pj0 = Rstat I02 = 22, 40W ± 2% 2 Quindi le perdite nel ferro e meccaniche si ricavano dalla relazione: Pf m = P0 − Pjo = (380, 80 ± 4, 89)W prova n° 1 2 3 4 5 6 7 8 9 10 V0 [V] I0 [A] P0 [W] cosϕ0 S0 [VA] Q0 [VAR] Pj0 [W] Pf m [W] 121,2 1,36 130,3 0,46 284,9 225,3 1,79 128,5 138,0 1,49 142,3 0,40 355,2 298,2 2,15 140,2 159,7 1,63 158,1 0,35 450,3 394,8 2,58 155,6 182,5 1,86 175,6 0,30 586,4 533,8 3,35 172,3 221,3 2,24 202,5 0,24 858,9 811,2 4,89 197,6 264,7 2,8 243,1 0,19 1264,9 1218,2 7,4 235,7 297,1 3,21 275,6 0,17 1651,8 1605,9 10,0 265,6 343,7 3,96 334,8 0,14 2358,0 2310,5 15,3 319,5 381,0 4,80 403,2 0,13 3164,9 3113,6 22,4 380,8 398,6 5,33 451,3 0,12 3678,4 3623,0 27,6 423,7 Table 16.4.: Caratteristiche a vuoto Nella tabella sono state riportate anche le potenze apparenti e reattive: So = 164 √ 3Vo Io di C.Bettini, A.Guiggiani e T.Lorenzetti 16.3 Circuto equivalente a vuoto Q0 = So sin ϕ0 16.2.1. Separazione delle perdite nel ferro e meccaniche Le perdite Pf m sono proporzionali a V 2 . Considerando costanti con la tensione le perdite meccaniche e dipendenti da V 2 le perdite nel ferro, si ottiene la seguente relazione: Pf m = Pf + Pm = KV 2 + Pm Tale relazione non è altro che un’equazione di una retta dove Pm è l’intercetta. Si ottiene quindi che le perdite meccaniche, identificabili nel punto in cui la retta interseca l’asse y, sono pari a 105 W (Figura 16.1). Figura 16.1.: Grafico delle perdite ferro-meccaniche in funzione di V 2 16.3. Circuto equivalente a vuoto Trascurando le perdite per effetto Joule nel rame statorico e osservando che per lo scorrimento s → 0 ⇒ R(s) → ∞, è possibile calcolare i valori del circuito equivalente di C.Bettini, A.Guiggiani e T.Lorenzetti 165 16.4 Grafici a vuoto, tenendo conto dei dati letti nella prova effettuata alla tensione nominale. Si ottiene così il circuito in Figura 16.2. Figura 16.2.: Circuito semplificato a vuoto Rm = 3Vo2 = 1080Ω ± 2, 1% Po Xm = 3V02 = 140Ω ± 2, 1% Q0 16.4. Grafici In Figura 16.3, Figura 16.4 e Figura 16.5 sono riportati i grafici ottenuti dalle misure effttuate nella prova a vuoto. 16.5. Conclusioni Osservando i grafici si nota che: • la potenza Pf m ha un andamento parabolico dovuto alle perdite nel ferro (vedi: Separazione delle perdite nel ferro e meccaniche); 166 di C.Bettini, A.Guiggiani e T.Lorenzetti 16.5 Conclusioni Figura 16.3.: Andamento delle potenze rispetto alla tensione nella prova a vuoto • il fattore di potenza diminuisce in quanto all’aumentare della tensione aumenta il comportamento reattivo del motore trifase (le perdite nel ferro sono preponderanti rispetto a quelle per effetto Joule). Confrontando le grandezze nominali del motore con quelle ricavate a vuoto alla tensione nominale si nota che: • Potenza nominale assorbita a vuoto: P0n = 403, 2W ; Pon % = 7, 3% • Corrente nominale assorbita a vuoto I0n = 4, 8A ; I0n % = 41% • Perdite meccaniche Pm = 105W ± 1, 5% • Perdite nel ferro nominali Pf n = 275, 8W di C.Bettini, A.Guiggiani e T.Lorenzetti 167 16.5 Conclusioni Ϭ͕ϱ Ϭ͕ϰϱ &ĂƚƚŽƌĞĚŝWŽƚĞŶnjĂ;ĐŽƐφ φͿ Ϭ͕ϰ Ϭ͕ϯϱ Ϭ͕ϯ Ϭ͕Ϯϱ Ϭ͕Ϯ Ϭ͕ϭϱ Ϭ͕ϭ Ϭ͕Ϭϱ Ϭ Ϭ ϱϬ ϭϬϬ ϭϱϬ ϮϬϬ ϮϱϬ ϯϬϬ ϯϱϬ ϰϬϬ ϰϱϬ sŽ Figura 16.4.: Andamento del fattore di potenza rispetto alla tensione nella prova a vuoto ϲ ϱ /Ž;Ϳ ϰ ϯ Ϯ ϭ Ϭ Ϭ ϱϬ ϭϬϬ ϭϱϬ ϮϬϬ ϮϱϬ ϯϬϬ ϯϱϬ ϰϬϬ ϰϱϬ sŽ Figura 16.5.: Andamento della corrente rispetto alla tensione nella prova a vuoto 168 di C.Bettini, A.Guiggiani e T.Lorenzetti 17. Prova in corto La prova in corto circuito è da effettuarsi a rotore bloccato alimentando il motore inizialmente con un valore di corrente pari al 110% circa della nominale ed andando a decrementarla fino ad arrivare a circa il 25% . La prova serve a misurare la Pcc che si può considerare con una certa approsimazione pari alle perdite nel rame. Essendo il rotore bloccato ne deriva che lo scorrimento (s) è unitario. 17.1. Tabelle delle misure prova n° 1 2 3 4 5 6 7 8 9 10 Vcc [V] Icc [A] Pcc [W] cosϕcc Qcc [VAR] Scc [VA] 102,80 13,40 998,60 0,419 1968,9 2385,9 98,00 12,53 883,20 0,415 1760,1 2126,9 91,20 11,67 716,30 0,389 1565,1 1843,4 83,43 10,01 497,60 0,344 1275,3 1446,5 72,36 8,20 362,10 0,352 899,9 1027,5 68,36 7,59 338,70 0,370 771,4 899,0 58,88 6,14 224,4 0,358 546,2 626,6 51,99 5,14 158,3 0,342 408,2 462,4 44,34 4,04 99,32 0,320 278,2 310,0 39,35 3,34 69,06 0,304 206,5 227,4 Table 17.1.: Dati ricavati dalle misure in corto circuito Nella tabella sono state riportate anche le potenze apparenti e reattive, ricavate come: Scc = √ 3Vcc Icc 169 17.2 Circuito equivalente in corto Qcc = Scc sin ϕcc Nella prova in corto circuito le perdite nel ferro possono essere trascurate, in quanto la tensione durante la prova assume valori molto inferiori rispetto a Vn , mentre le perdite meccaniche sono nulle. Quindi: 2 Pcc = Pj1 + Pj2 = 3Reqcc Icc dove Pj1 e Pj2 sono le perdite per effetto Joule nel rame rispettivamente dello statore e del rotore ed Reqcc è la resistenza equivalente riportata allo statore. 17.2. Circuito equivalente in corto Nella prova in corto lo scorrimento è s = 1, quindi la resistenza fittizia R(s) è nulla. Con una certa approssimazione si può trascurare l’ammettenza Zm dato il ridotto valore di tensione della prova. Quindi si ottiene il seguente circuito semplificato in Figura 17.1 dove: Vccn Zeqcc = √ 3In Reqcc = Zeqcc cos ϕcc Xeqcc = Zeqcc sin ϕcc 170 di C.Bettini, A.Guiggiani e T.Lorenzetti 17.3 Parametri riportati alla temperatura convenzionale sono le grandezze equivalenti riportate allo statore essendo Vccn la tensione di corto circuito riferita alla corrente nominale. Figura 17.1.: Circuito equivalente in corto circuito (semplificato) 17.3. Parametri riportati alla temperatura convenzionale I parametri caratteristici della prova in corto circuito vanno riportati alla temperatura convenzionale di Tr = 75°; è importante notare che la reattanza non dipende dalla temperatura, quindi: T Reqcc = K T Reqcc T = Xeqcc Xeqcc T Zeqcc = q T )2 + (X T )2 (Reqcc eqcc di C.Bettini, A.Guiggiani e T.Lorenzetti 171 17.4 Grafici dove K T si trova con la formula ??eq:riporto a pagina 160. Il fattore di potenza alla temperatura di riferimento è pari a: cos ϕTcc = T Reqcc T Zeqcc (17.1) mentre la potenza persa in corto circuito alla temperatura di riferimento è: PccT = K T 3Reqcc In2 = K T Pcc = 859W ± 1, 9% Tp = 23, 5° [Ω] Accuracy Tr = 75° T Reqcc 1,85 ±4, 4% Reqcc T Xeqcc 3,76 ±4, 4% Xeqcc T Zeqcc 4,51 ±1, 1% Zeqcc Tabella 17.3.: Parametri in corto [Ω] Accuracy 2,21 ±4, 9% 3,76 ±4, 4% 4,36 ±4, 5% circuito 17.4. Grafici In Figura 17.2, Figura 17.3 e Figura 17.4 sono riportati i grafici delle misure effettuate nella prova in corto circuito. 17.5. Conclusioni Andando ad osservare i grafici, si nota che: • La Vcc ha un andamento lineare rispetto alla corrente in quanto il comportamento è legato alla caduta di tensione sull’impedenza equivalente secondo la √ relazione lineare Vcc = 3Zeqcc Icc ; 172 di C.Bettini, A.Guiggiani e T.Lorenzetti 17.5 Conclusioni T Figura 17.2.: Andamento della tensione rispetto alla corrente in corto circuito • Le perdite nel rame Pcc hanno un andamento parabolico essendo dipendenti dal quadrato della corrente; • Il fattore di potenza cos ϕcc si mantiene quasi costante in quanto non è dipen eqcc dente dalla corrente cos ϕcc = R . Zeqcc Confrontando le grandezze nominali del motore con quelle ricavate in corto circuito (alla corrente nominale) si nota che: • Potenza nominale assorbita in corto circuito: Pcc = 716, 3W ; Pccn % = 13% • Corrente nominale assorbita in corto circuito Vccn = 91, 2V ; di C.Bettini, A.Guiggiani e T.Lorenzetti Vccn % = 24% 173 17.5 Conclusioni Figura 17.3.: Andamento delle perdite nel rame rispetto alla corrente in corto circuito Figura 17.4.: Andamento del fattore di potenza rispetto alla corrente in corto circuito 174 di C.Bettini, A.Guiggiani e T.Lorenzetti 18. Diagramma circolare Dalle prove effettuate possiamo ricavare i dati necessari per la costruzione del diagramma circolare. Il diagramma circolare è un diagramma polare che mostra il variare della corrente assorbita in funzione dello scorrimento supponendo la tensione di alimentazione costante. Inoltre è possibile anche determinare le potenze, assorbite e rese, con le relative perdite e la coppia in funzione dello scorrimento. Per ottenere questi paramentri è necessario calcorare le intersezioni fra le relative rette: • L’asse delle ordinate rappresenta la potenza assorbita dal motore; • La retta “a” rappresenta le potenze rese all’asse del motore; • La retta “b” rappresenta le coppie trasmesse. Prendendo come esempio un punto di lavoro con scorrimento s = 0, 033, possiamo ricavarci i parametri valutando i seguenti segmenti: AB Potenza assorbita BD Potenza persa nello statore DE Potenza dissipata nel rotore AD Potenza trasmessa al rotore AE Potenza resa Quindi il rendimento del motore nel punto di riferimento è dato da; η= AE AB = Pr Pa 175 Diagramma circolare Dal grafico (Figura 18.1 e Figura 18.2) si ricavano le seguenti misure: AB = 8, 7 AE = 7, 27 AD = 7, 72 BD = 0, 98 DE = 0, 45 n2 = 1500 − (1500 · 0, 033) = 1450rpm Per ricavare i valori delle potenze è necessario moltiplicare i segmenti per il fattore √ di scala 3Vn . Si ottiene: AB = 5, 726kW =Potenza assorbita AE = 4, 785kW =Potenza resa AD = 5, 081kW =Potenza trasmessa al rotore BD = 645W =Potenza dissipata nello statore DE = 296W =Potenza dissipata nel rotore Ne deriva un rendimento pari a: η= 176 AE AB ' 83, 57% di C.Bettini, A.Guiggiani e T.Lorenzetti Diagramma circolare di C.Bettini, A.Guiggiani e T.Lorenzetti Figura 18.1.: Diagramma circolare con s = 0, 033 177 Diagramma circolare Figura 18.2.: Particolare del diagramma circolare del caso in esame (s = 0, 033) 178 di C.Bettini, A.Guiggiani e T.Lorenzetti Parte V. Laboratorio di automatica I 179 19. Controllo di un levitatore magnetico con PID e tecnica di taratura automatica IFT Nella seguente relazione è trattata l’applicazione di un controllo iterativo non model based ad un levitatore magnetico. Il sistema in oggetto è descritto all’interno del Capitolo 2 ed è composto da un elettromagnete che agisce su una sfera metallica inducendo una forza su di essa. Il nostro scopo è quello di riuscire a controllare il suddetto sistema inseguendo varie tipologie di riferimento. Il sistema a disposizione è comandato da una macchina con scheda di acquisizione ed elaborazione Real Time, quindi è stato necessario lavorare in tempo discreto. Per controllare il sistema è stato utilizzato un controllore di tipo PID il quale è stato tarato tramite l’algoritmo IFT (Iterative Feedback Tuning). L’IFT è un metodo di taratura iterativo che cerca di migliorare le prestazioni del sistema andando a minimizzare un determinato funzionale di costo. Per calcolare i parametri ottimi del controllore PID, l’IFT presuppone l’esistenza di due esperimenti durante i quali vengono acquisiti i segnali di uscita, controllo ed errore che verranno poi elaborati dall’algoritmo. La problematica affrontata è stata quindi quella di riuscire a temporizzare lo schema Simulink Real-Time in modo tale da poter eseguire gli esperimenti non solo in modo iterativo ma anche in modo consequenziale; infatti lo schema Simulink deve implementrare due esperimenti sul sistema in anello chiuso: il primo con in ingresso il segnale di riferimento ed il secondo con ingresso l’errore del primo esperimento. 181 Controllo di un levitatore magnetico con PID e tecnica di taratura automatica IFT Inoltre ad ogni iterazione è necessario effettuare un aggiornamento dei parametri (opportunamente inizializzati) secondo la relazione: pi+1 = pi − ∆p in cui l’incremento ∆p è calcolato in base alla minimizzazione di una determinata funzione di costo nel nostro dipendente dai segnali di errore e di uscita del sistema in anello chiuso. Al termine dell’analisi sono stati ricavati dei PID ottimi per i diversi segnali di riferimento presi in esame. 182 di C.Bettini, A.Guiggiani e T.Lorenzetti 20. Modello In questo capitolo si cerca di definire ed identificare il modello fisico del processo partendo da un’analisi diretta delle forze che agiscono sul sistema. Tutto ciò è atto ad ottenere un modello matematico sufficientemente coerente con quello reale in modo da poter valutare in fase di simulazione le prestazioni del controllore PID che andremo a migliorare. Il modello matematico complessivo che descrive il processo è composto dall’unione dei comportamenti dei singoli dispositivi del sistema: un amplificatore di transconduttanza, un elettromagnete ed un sensore di posizione. 20.1. Amplificatore di transconduttanza L’amplificatore di transconduttanza è un attuatore, ovvero un elemento che fornisce un segnale di potenza (nel caso specifico la corrente di eccitazione nelle spire dell’elettromagnete) proporzionale al segnale di controllo in tensione inviato in ingresso. Ipotizzeremo che nel nostro caso esso abbia (indipendentemente dalle condizioni di lavoro) banda infinita; ciò significa che l’amplificatore è considerato a tutti gli effetti ideale e caratterizzato da un tempo di risposta nullo. Inoltre, visto che la forza magnetica è solo attrattiva (dato che il magnete è posizionato sopra alla sfera), nella caratteristica statica dell’amplificatore troveremo un valore di corrente nullo per tensioni di ingresso negative mentre - per ragioni di sicurezza - è presente un circuito di protezione che azzera la corrente in ingresso all’elettromagnete quando questa supera la soglia di 3 A. 183 20.2 Elettromagnete Questo permette di descrivere il comportamento dell’amplificatore attraverso l’espressione: Im = Ka Vin + i0 ∀ 0 ≤ Im ≤ 3 (20.1) in cui si ha: Im =corrente in uscita dell’amplificatore Ka = 0.4087[S] =guadagno dell’amplificatore Vin =tensione di ingresso i0 = 0.0944[A] =corrente di offset corrispondente a Vin = 0 20.2. Elettromagnete L’elettromagnete è costituito da un nucleo sagomato ad E sul cui elemento centrale sono avvolte le spire che percorse da corrente generano il campo magnetico indotto. Dalle formule fondamentali dell’elettromagnetismo si ottiene l’equazione dell’energia accumulata nel traferro: ˆ ˆ W = H ! B dH dV V 0 (20.2) mentre la forza di attrazione magnetica è pari a: F =− dW dz (20.3) dove: V =volume al traferro 184 di C.Bettini, A.Guiggiani e T.Lorenzetti 20.2 Elettromagnete H =intensità del campo magnetico B =induzione magnetica z =coordinata spaziale in cui si calcola la forza quindi, logicamente, a parità di corrente nell’elettromagnete l’intensità del campo magnetico H (e di conseguenza della forza F ) varia in funzione della posizione z della sfera. Per sistemi a geometria molto semplice, ipotizzando per semplificare che la permeabilità del nucleo sia infinitamente più grande di quella del vuoto, la forza con cui viene attratta la sfera è derivabile dalle equazioni: N I = 2lH B = µ0 H W = µ0 H 2 Sl = F = µ0 SN 2 1 µ0 SN 2 I 2 4l I2 I2 = K (2l)2 (2l)2 dove: µ0 =permeabilità magnetica del vuoto N =numero di spire dell’avvolgimento I =corrente nelle spire dell’avvolgimento S =sezione del traferro di C.Bettini, A.Guiggiani e T.Lorenzetti 185 20.3 Trasduttore di posizione 2l = lunghezza totale del traferro Dato che la permeabilità del nucleo viene considerata infinitamente più grande di quella del vuoto, allora si può ipotizzare che la forza magnetomotrice NI sia concentrata completamente nel traferro e che non sia presente flusso disperso (ovvero tutte le linee di flusso passano per il traferro). L’aver supposto di considerare elettromagneti con geometria semplice non toglie validità alle considerazioni già fatte; risulta infatti ragionevole ipotizzare che tra la forza magnetica Fm , la corrente Im e la distanza z esista una relazione esprimibile come: 2 Im (t) Fm (t) = Km 2 z (t) (20.4) con: Km = 1.6284 · 10−4 h N m2 A2 i = costante elettromagnetica Tenendo conto che sulla sfera agiscono sia la forza magnetica Fm sia la forza peso e trascurando i fenomeni di attrito viscoso, si può facilmente trovare l’equazione che descrive la dinamica del moto: ms z̈(t) = ms g − Fm (t) = ms g − Km 2 (t) Im 2 z (t) (20.5) con ms pari alla massa della sfera e g all’accelerazione gravitazionale. 20.3. Trasduttore di posizione Il trasduttore di posizione è un dispositivo di tipo ottico realizzato in configurazione differenziale, cioè in modo tale da discriminare i movimenti verticali da quelli orizzontali della pallina. 186 di C.Bettini, A.Guiggiani e T.Lorenzetti 20.4 Modello matematico del sistema Il segnale luminoso emesso dai fotodiodi non è a intensità costante ma costituito da una portante a 50 kHz ed ampiezza costante. I segnali in uscita dal ricevitore, anch’essi sinusoidali a 50 kHz, sono modulati in ampiezza dalla quantità d’ombra proiettata dalla sfera e successivamente demodulati prima di essere opportunamente elaborati al fine di produrre il segnale di uscita Vout . Tale realizzazione è stata scelta al fine di rigettare i disturbi luminosi provenienti dall’ambiente circostante. Nonostante la sua complessità, il trasduttore può essere considerato statico in quanto la sua velocità di risposta è di alcuni ordini di grandezza superiore a quella delle restanti parti del sistema con un comportamento non lineare man mano che ci avviciniamo agli estremi del suo campo di azione. Nell’intervallo di linearità questo dispositivo può essere adeguatamente descritto dalla equazione: Vout = Kt z + V0 (20.6) dove: Kt = 529.9527 h i V m =guadagno di transconduttanza V0 = −12.0715[V ] =tensione di offset 20.4. Modello matematico del sistema Dalle equazioni (Equazione 21.1), (Equazione 21.3) e (Equazione 21.2)si ricava il modello matematico del sistema: (t)+i0 ) ms z̈(t) = ms g − Fm (t) = ms g − Km (Ka Vin z 2 (t) Vout = Kt z + V0 2 (20.7) rappresentato dal diagramma a blocchi in Figura 34.1. di C.Bettini, A.Guiggiani e T.Lorenzetti 187 20.5 Tempo di campionamento Figura 20.1.: Schema Simulink del processo tempo-discreto Tutti i valori delle costanti caratterizzanti il modello sono state ricavate in modo empirico da uno studio precedente del sistema fisico. 20.5. Tempo di campionamento Occorre adesso scegliere il tempo di campionamento Ts da utilizzare per implementare il controllore in tempo discreto. La scelta di Ts non è particolarmente problematica in quanto la scheda di acquisizione usata consente la possibilità di lavorare con frequenze di campionamento molto elevate rispetto a quella che è la velocita di risposta dell’impianto con cui si ha a che fare. La scelta ritenuta più opportuna è stata quella di imporre Ts = 0.001[s] al fine di apprezzare le dinamiche più veloci del processo senza alcun problema di ricostruibilità dei segnali campionati risultando molto affini al corrispondente segnale in tempo continuo. 188 di C.Bettini, A.Guiggiani e T.Lorenzetti 21. Cenni sui regolatori PID Iniziamo con il descrivere il funzionamento dei regolatori lineari con più ampia applicazione in ambito industriale: i PID (Proporzionale-Integrativo-Derivativo). Un così ampio successo è dovuto a ragioni molteplici: • Il loro impiego consente di controllare una ampia gamma di processi; • Semplice costruzione: questa peculiarità permette ai PID di essere realizzati con le tecnologie più varie (meccanica, idraulica, elettronica digitale o analogica, ecc); • Negli anni sono state sviluppate diverse tecniche di taratura automatica che permettono di controllare il sistema di interesse senza conoscerne l’esatto modello matematico. 21.1. Il modello Come si intuisce dal nome, la variabile di controllo u generata da regolatore di tipo PID sarà composta da tre contributi: Proporzionale: è un contributo proporzionale rispetto all’errore (e) tra il segnale di riferimento in ingresso r e l’uscita y. Tale contributo non influisce sulla stabilità, ma migliora le prestazioni statiche del sistema da controllare. Integrativo: è un contributo proporzionale all’integrale di e. Esso può essere interpretato come una azione ritardatrice sul sistema, che va a soddisfare i requisiti sull’errore a transitorio esaurito in risposta ad un riferimento a gradino. Derivativo: è un contributo proporzionale alla derivata di e. Esso viene interpretato come una azione anticipatrice, utile in quei casi in cui sia necessario ottenere 189 21.1 Il modello più banda passante possibile, andando quindi a migliorare le prestazioni del transitorio del sistema. La legge di controllo che ne deriva è la seguente: ˆ u(t) = Kp e(t) + Ki e(τ )dτ + Kd de(t) dt (21.1) Dove Kp , Ki e Kd sono rispettivamente i coefficienti dell’azione proporzionale, integrale e derivativa. Si dimostra facilmente che per i PID, almeno nella loro forma ideale, sono verificate tutte quelle ipotesi che ci garantiscono la correttezza dell’uso della trasformata di Laplace. Da ciò ne deriviamo quindi la forma di più facile utilizzo: CP ID (s) = KP + KI + Kd s s (21.2) La forma che useremo noi, forse anche la più usata è la seguente : CP ID (s) = KP 1 1+ + Td s Ti s (21.3) dove Ti = Kp /Ki è il tempo integrale e Td = Kd /Kp è il tempo derivativo. La funzione di trasferimento del PID nella sua forma ideale è un sistema improprio, in quanto avente un numero di zeri maggiore ai poli a causa del termine derivativo. Per questo motivo, nella pratica viene utilizzata la seguente funzione di trasferimento: CP ID (s) = KP Td 1 1+ + Ti s 1 + TNd s ! (21.4) Dove la costante N viene scelta propriamente in modo che il polo aggiunto s = − TNd sia fuori dalla banda di controllo. 190 di C.Bettini, A.Guiggiani e T.Lorenzetti 21.2 Realizzazione dei regolatori PID Le tipiche rappresentazioni di Bode dei PID sono le seguenti (sia reale che ideale): Figura 21.1.: Diagrammi di Bode asintotici ed effettivi di un PID ideale e reale 21.2. Realizzazione dei regolatori PID Lo schema di partenza di un controllore PID può essere visto come un parallelo di tre blocchi corrispondenti alle tre azioni (proporzionale, derivativo, integrale) del tipo rappresentato in Figura 21.2. Figura 21.2.: Schema di controllo con regolatore PID a derivazione dell’errore Sì può notare però che applicando direttamente l’azione derivativa sull’errore può accadere che, mettendo particolari segnali in ingresso (per esempio un gradino) al sistema, si ottenga una azione di controllo u di natura impulsiva. Questa brusca reazione va contro il requisito di moderazione del controllo e può provocare la saturazione dell’attuatore e l’allontanamento del sistema dalla linearità. Quindi viene di C.Bettini, A.Guiggiani e T.Lorenzetti 191 21.2 Realizzazione dei regolatori PID usato uno schema alternativo riportato in Figura 21.3, dove l’azione derivativa viene applicata sull’uscita e non sull’errore, poiché è giusto ritenere che il processo abbia una caratteristica di un filtro passa-basso assicurando un’uscita senza brusche variazioni (e quindi la sua derivata non ha una caratteristica impulsiva). Figura 21.3.: Schema di controllo con regolatore PID a derivazione dell’uscita Per quanto riguarda il parametro N il progettista deve cercare di aumentare N per mettere il polo aggiunto s = − TNd fuori banda di controllo avvicinandosi così al comportamento ideale del PID. D’altra parte, deve scegliere anche un valore di N sufficientemente basso per cui i disturbi in alta frequenza non vengano eccessivamente amplificati. Valori tipici di N vanno da 5 a 20. 192 di C.Bettini, A.Guiggiani e T.Lorenzetti 22. Iterative Feedback Tuning Generalmente, soluzioni esplicite di obiettivi di controllo richiedono una completa conoscenza del sistema dalle caratteristiche frequenziali e/o probabilistiche, dei disturbi e la possibilità di lavorare con controllori arbitrariamente complessi. Nella pratica spesso il modello del processo e le perturbazioni non sono noti e si preferisce lavorare con controllori di prestabilita complessità, come i PID. Per raggiungere le migliori prestazioni, ci si affida a metodi iterativi basati sulla minimizzazione di una funzione di costo. Il problema principale quindi con queste tecniche è relazionare il gradiente della funzione di costo ai parametri del controllore, poiché se il modello e i disturbi non sono conosciuti non è chiaro come il gradiente possa essere calcolato. In primo luogo sì cercò di determinare metodi di identificazione ad anello chiuso che garantissero buone prestazioni. Questo fu possibile solo su modelli complessi (fullorder models) utilizzando controllori complessi, situazione ideale però irrealistica. Usando questi modelli è stato così possibile dimostrare che attraverso questo tipo di identificazione si ottiene le prestazioni ottime sul sistema. Nel caso di controllori semplici (low-order controllers) sono stati riportati numerosi successi, seguendo gli schemi sopracitati basati sull’identificazione del modello. Si notò però che in alcuni casi non si ottiene la convergenza, per cui si pensò che la strada migliore da seguire fosse la determinazione diretta dei parametri del controllore (piuttosto che del modello del sistema). Comunque sia, come è stato detto precedentemente, le prime prove in questa direzione si sono scontrate contro la difficoltà di calcolare il gradiente del costo rispetto ai parametri del controllore. Dagli studi si è concluso che il metodo migliore per stimare il gradiente deriva dalla rielaborazione dei segnali ottenuti da esperimenti ad anello chiuso utilizzando il sistema dato e il controllore. Il numero degli esperimenti dipende strettamente dalla 193 Iterative Feedback Tuning complessità del controllore. Per controllori a due gradi di libertà, sono sufficienti tre esperimenti dove il primo e il terzo servono a registrare semplicemente l’evoluzione dei segnali; mentre il secondo (che è il vero e proprio esperimento) utilizza come riferimento di ingresso i dati raccolti dal primo esperimento. Proprio su questi esperimenti si basa l’Iterative Feedback Tuning (IFT), presentato la prima volta nel 1994, dove è parso chiaro, data la semplicità dello stesso, il vasto potenziale applicativo. In particolare l’attrattiva di tale metodo in ambito ingegneristico è dovuta alla possibilità di ottimizzare i parametri del controllore osservando il comportamento del sistema ad anello chiuso. Inoltre con questo schema, la regolazione dei parametri con reiezione dei disturbi è guidata dai disturbi stessi, rendendo molto semplice quindi l’operare in presenza di processi stocastici. Consideriamo un sistema sconosciuto descritto da un modello a tempo continuo: Figura 22.1.: Sistema con un controllore a due gradi di libertà y = Gu + v (22.1) dove G è un sistema LTI, v è un disturbo (processo) aleatorio, y l’uscita del sistema e u l’ingresso del sistema LTI. Inoltre assumiamo che il sistema sia controllato da un controllore a due gradi di libertà : u = Cr (r)r − Cy (r)y 194 di C.Bettini, A.Guiggiani e T.Lorenzetti (22.2) Iterative Feedback Tuning dove Cr (r) e Cy (r) sono funzioni di trasferimento LTI parametrizzate da un vettore di parametri r Rn e r è un segnale di riferimento deterministico indipendente dai disturbi. Sia y d l’uscita desiderata ottenuta applicando il segnale di riferimento r al sistema ad anello chiuso che abbia una funzione di trasferimento desiderata Td ovvero: y d = Td r (22.3) Quindi l’errore è definito come: e(r) = y(r)–y d = 1 Cr (r)G r − yd + v 1 + Cy (r)G 1 + Cy (r)G (22.4) oppure se il modello di riferimento Td è definito si può scrivere anche: Cr (r)G 1 –Td r + v 1 + Cy (r)G 1 + Cy (r)G ! e(r) = (22.5) L’errore è dato da due contributi: • dal non corretto inseguimento del segnale di riferimento r; • dalla componente di disturbo. L’Iterative Feedback Tuning (IFT) affronta il problema del controllo ottimo come un problema di minimizzazione del gradiente di una funzione di costo, nel nostro caso dipendente dall’errore e(r) e dal controllo u(r). Consideriamo la seguente forma quadratica: J(r) = E (ˆ T 0 2 2 ) e (t) + lu (t) dt di C.Bettini, A.Guiggiani e T.Lorenzetti (22.6) 195 22.1 Criteri di minimizzazione I parametri ottimi r sono quelli che minimizzano tale funzione di costo. r∗ = arg min J(r) (22.7) ρ Il primo termine nella funzione di costo J(ρ) è proporzionale all’errore del sistema, mentre il secondo, pesato su l, è proporzionale ai controlli. 22.1. Criteri di minimizzazione Vediamo ora la minimizzazione della funzione di costo J(r) (Equazione 22.6) rispetto al vettore dei parametri r di un controllore con una struttura già preimpostata. Appare evidente a partire dalla (Equazione 22.4), che J(r) dipende in un modo abbastanza complicato da r, dal sistema sconosciuto G e dallo spettro sconosciuto di v. Per ottenere una soluzione di minimo di J(r) dobbiamo trovare una soluzione alla seguente espressione: ∂J (r) = E 2 ∂r ˆT ! e(t) 0 ∂e ∂u + lu(t) dt = 0 ∂r ∂r (22.8) Supponendo di poter calcolare il gradiente ∂J (r), sarebbe possibile risolvere la ∂r (Equazione 22.8) attraverso un algoritmo iterativo in cui il vettore dei parametri è aggiornato ricorsivamente secondo la relazione: ri+1 = ri − gi Ri−1 ∂J (r ) ∂r i (22.9) dove Ri è una matrice definita positiva, tipicamente una approssimazione GaussNewton della matrice Hessiana di J; mentre gi è uno scalare reale e positivo che 196 di C.Bettini, A.Guiggiani e T.Lorenzetti 22.1 Criteri di minimizzazione determina il passo di iterazione. La sequenza gi deve rispettare dei vincoli perché l’algoritmo converga ad un minimo locale della funzione di costo J. (ri ) è intrattabile a meno che non si ricorra ad un Il problema del calcolo di ∂J ∂r algoritmo di approssimazione stocastica come suggerito da Robbins e Monro[?]. Il modo per risolvere questo problema è il seguente: 1. Memorizzare e(ri ) e u(ri ); 2. Stimare il gradiente ∂e (r ) ∂r i e ∂u (r ); ∂r i 3. Determinare le stime dei prodotti e(ri ) ∂∂er (ri ) e u(ri ) ∂u (r ). ∂r i La stima dei due gradienti è sempre stato l’ostacolo per risolvere il problema. Tuttavia, come si vedrà in seguito queste quantità possono essere stimate facendo degli esperimenti sul sistema ad anello chiuso con i controllori Cr (r) e Cy (r). 22.1.1. Calcolo dei gradienti Dalla (Equazione 22.4) nella pagina 195 è chiaro che e(r) è ottenuto dalla differenza tra la risposta con il controllore C(r) = {Cr (r), Cy (r)} e quella desiderata. Si può (r). Quindi: notare anche che ∂∂er (r) = ∂y ∂r G ∂Cr Cr (r)G2 ∂Cy G ∂Cy ∂e (r) = (r)r − (r)r − (r)v 2 2 ∂r 1 + Cy (r)G ∂ r 1 + Cy (r)G) ∂ r (1 + Cy (r)G) ∂ r (22.10) ∂y 1 ∂Cr 1 ∂Cy (r) = (r)T r– (r)(T 2 r + T Sv) ∂r Cr (r) ∂ r (Cr (r) ∂ r (22.11) In questa espressione le quantità 1/Cr (r), ∂Cr /∂ r e∂Cy /∂ r sono funzioni conosciute, mentre T e S sono funzioni di un sistema incognito e quindi non sono calcolabili. Per calcolare il segnale ∂y/∂ r si ricorre a degli esperimenti sul sistema ad anello chiuso. di C.Bettini, A.Guiggiani e T.Lorenzetti 197 22.1 Criteri di minimizzazione Ora osserviamo che l’ultimo termine può essere riscritto come: T 2 (r)r + T (r)S(r)v = T (r)y (22.12) Quindi la (Equazione 22.11) si può scrivere così: ∂y 1 (r) = ∂r Cr (r) " ! # ∂Cr ∂Cy ∂Cy (r) − (r) T (r)r + (r)T (r)(r − y) ∂r ∂r ∂r (22.13) L’ultimo termine può essere ottenuto sottraendo il segnale di uscita dal segnale di riferimento ottenuto dall’esperimento, usando questo segnale d’errore come segnale di riferimento nell’esperimento successivo. Questa osservazione ci induce a suggerire la seguente procedura. In ogni iterazione i-esima dell’algoritmo faremo due esperimenti della stessa durata. Prendiamo un segnale di riferimento rij (con j =1,2 che indica il numero dell’esperimento e i la i-esima iterazione) e il corrispondente segnale di uscita y j (r). Avremo che : I esperimento ri1 = r y 1 (r) = T (r)r + S(r)vi1 II esperimento ri2 = r − y 1 (r) y 2 (r) = T (r)(r − y 1 (r)) + S(r)vi2 I disturbi possono essere considerati mutuamente indipendenti, venendo da differenti esperimenti. Quindi, sostituendo possiamo ora scrivere una stima del gradiente: ∂y 1 est[ (r)] = ∂r Cr (r) 198 " ! # ∂Cr ∂Cy ∂Cy (r) − (r) y 1 (r) + (r)y 2 (r) ∂r ∂r ∂r di C.Bettini, A.Guiggiani e T.Lorenzetti (22.14) 22.1 Criteri di minimizzazione Ora vediamo come si arriva al calcolo della stima del gradiente dei controlli : u(r) = Cr (r) Cy (r) r+ v 1 + Cy (r)G 1 + Cy (r)G (22.15) ∂u ∂Cr ∂S ∂Cy ∂S (r) = (r)S(r)r − (r)Cr (r)r − (r)S(r)v − (r)Cy (r)v (22.16) ∂r ∂r ∂r ∂r ∂r 1 ∂Cy ∂S (r) = S(r)T (r) (r) ∂r Cr (r) ∂r Sostituendo : " # ∂Cr ∂Cy ∂u (r) = S(r) (r)r − (r)(T (r)r − S(r)v) ∂r ∂r ∂r (22.17) T (r)r + S(r)v = y (22.18) " ∂u ∂Cr ∂Cy (r) = S(r) (r)r − (r)y ∂r ∂r ∂r Aggiungendo e togliendo " # (22.19) ∂Cy (r)r: ∂r ∂u ∂Cr ∂Cy ∂Cy (r) = S(r) (r)r − (r)r + (r)(r − y) ∂r ∂r ∂r ∂r # di C.Bettini, A.Guiggiani e T.Lorenzetti (22.20) 199 22.2 Stima del gradiente Come in precedenza in questa espressione le quantità ∂Cr/∂ r(r) e ∂Cy/∂ r(r) sono funzioni conosciute, mentre S è una funzione incognita, quindi per calcolare il segnale ∂u/∂ r(r) si ricorre ai due esperimenti sul sistema ad anello chiuso. I esperimento ri1 = r u1 (r) = S(r) [Cr (r)r − Cy (r)vi1 ] II esperimento ri2 = r − y 1 (r) u2 (r) = S(r) [Cr (r)ri2 − C(r)vi2 ] Sostituendo si può ottenere una stima del gradiente : ∂u 1 est (r) = ∂r Cr (r) " # " ! # ∂Cr ∂Cy ∂Cy (r) − (r) u1 (r) + (r)u2 (r) ∂r ∂r ∂r (22.21) 22.2. Stima del gradiente Partendo dalla (Equazione 22.8) nella pagina 196, possiamo stimare la derivata del gradiente della funzione di costo, ottenendo: # " ∂J est (r) = 2 ∂r ˆT ( 0 " # " ∂y ∂u e(t)est (r) + lu(t)est (r) ∂r ∂r #) dt (22.22) Perché l’approssimazione stocastica dell’algoritmo funzioni, è necessario che la stima del gradiente non sia polarizzata, e cioè: ( " #) ∂J E est (r) ∂r = ∂J (r) ∂r (22.23) Si verifica che tale condizione è rispettata, poiché gli errori usati per la stima del gradiente sono incorrelati (venendo da esperimenti, come detto prima, con disturbi mutuamente indipendenti). 200 di C.Bettini, A.Guiggiani e T.Lorenzetti 22.3 Convergenza 22.3. Convergenza In questo paragrafo discuteremo le condizioni di convergenza dell’algoritmo. Prendiamo in considerazione un sottoinsieme di Rn convesso e compatto D. Le condizioni per cui l’IFT converga ad un punto stazionario sono: 1. In ogni esperimento il segnale di disturbo v è un processo casuale e limitato per tutta la durata dell’esperimento |v| ≤ C ∀t. La costante C e le statistiche del secondo ordine del segnale v sono uguali in tutti gli esperimenti, ma i segnali dai differenti esperimenti sono mutuamente indipendenti. 2. Deve esistere un intorno ō di D tale che il controllore C(r) sia due volte differenziabile rispetto ai parametri r in ō. 3. Tutti gli elementi delle funzioni di trasferimento Cr (r) e Cy (r) e le loro derivate parziali hanno i poli e gli zeri con parte reale uniformemente minore di zero in D. 4. Il sistema LTI in anello chiuso è stabile e ha i poli con parte reale uniformemente minore di zero in D. 5. Gli elementi della sequenza gi devono essere positivi e devono rispettare i P∞ 2 P seguenti vincoli: ∞ i=1 gi < ∞. i=1 gi = ∞ ; Teorema: Supponiamo le condizioni sopra elencate vere. Supponiamo Ri sia una matrice generata da un esperimento alla iterazione i che soddisfa 1/dI ≥ Ri ≥ dI con d > 0. Allora per i → ∞, ρi converge con probabilità 1 a un punto dell’insieme: ( ) ∂J (r) = 0 Dc , r : ∂r (22.24) 22.4. Implementazione Per implementare il controllo abbiamo deciso di usare il regolatore PID. Tale controllore come è noto dipende da tre parametri: Kp , Ti e Td dove la relativa funzione di trasferimento è : di C.Bettini, A.Guiggiani e T.Lorenzetti 201 22.4 Implementazione CP ID (s) = KP 1 sTd 1+ + Ti s 1 + TNd s ! La struttura del sistema in retroazione da noi usata è la seguente: Figura 22.2.: Schema del sistema con il controllore PID con derivazione sull’uscita dove l’azione derivativa è posta sull’uscita per avere una moderazione del controllo rispetto a certi segnali di riferimento impulsivi, come il gradino. Tuttavia, uno schema di questo tipo per i nostri scopi è poco maneggevole in quanto non corrisponde al modello considerato nella trattazione teorica. Considerando quindi che il controllore trattato è a due gradi di libertà, possiamo scrivere: u = (r − y)CP I − yCD = rCP I − yCP I − yCD (22.25) u = rCP I − y(CP I + CD ) = rCP I − yCP ID (22.26) Si può quindi vedere che lo schema in Figura 22.2 è equivalente a quello in Figura 22.3. Questo schema è di più semplice implementazione e ci permette di usare direttamente le formule già viste in precedenza. Inoltre il controllore sul riferimento non 202 di C.Bettini, A.Guiggiani e T.Lorenzetti 22.4 Implementazione Figura 22.3.: Schema con controllore a due gradi di libertà contiene l’azione derivativa, così in presenza del gradino, come segnale di riferimento, il segnale di controllo u rispetta i requisiti di moderazione del controllo. Le funzioni di trasferimento dei due controllori sono le seguenti: CP I (s) = KP 1 1+ Ti s Td 1 + 1+ Ti s 1 + TNd s CP ID (s) = KP (22.27) ! (22.28) Dove Kp ,Ti e Td sono i parametri di controllo quindi il vettore r è costituito da K p ρ = Ti Td Come già accennato nel capitolo precedente, la stima del gradiente è il punto focale dell’algoritmo, ed è definita come: " # ∂J est (r) = 2 ∂r # ˆT " ∂y ∂u −e(t)est[ (r)] + lu(t)est[ (r)] dt ∂r ∂r (22.29) 0 Come si può notare, la stima del gradiente dipende dalle stime delle derivate parziali, sia di y e di u, rispetto ai parametri. Nel nostro caso abbiamo un controllore a tre di C.Bettini, A.Guiggiani e T.Lorenzetti 203 22.4 Implementazione parametri, per cui la stima del gradiente diventa: i ´T h ∂y ∂u ∂J 2 −e(t)est[ ( r )] + l u(t)est[ ( r )] dt est[ ( r )] ∂KP ∂KP ∂KP i ´0 T h ∂J ∂y ∂u ∂J (r) = est[ ∂Ti (r)] = 2 0 −e(t)est[ ∂Ti (r)] + lu(t)est[ ∂Ti (r)] dt est i ∂r ´T h ∂y ∂J ∂u est[ ∂T ( r)] ( r )] + l u(t)est[ ( r )] dt 2 −e(t)est[ ∂Td ∂Td 0 d # " Da questa formula appare chiara la necessità di calcolare quindi tutte le stime delle derivate parziali. Quindi sostituendo otteniamo le stime dei gradienti sull’uscita e sul controllo: " est # ∂y (r) = ∂r ∂y est[ ∂K (r)] P ∂y est[ ∂Ti (r)] ∂y est[ ∂T (r)] d − − = ∂CP I (r) ∂KP ∂CP I (r) ∂Ti ∂CP I (r) ∂Td ∂u PI (r) − est[ ∂K (r)] ∂C ∂KP P ∂u ∂u PI est (r) = est[ ∂T = ∂C (r) − (r)] i ∂Td ∂r ∂u ∂CP I est[ ∂Td (r)] (r) − ∂Td " # − ∂CP ID (r) y 1 (r) ∂KP ∂CP ID ( r ) y 1 (r) ∂K ∂CP ID (r) y 1 (r) ∂K ∂CP ID (r) u1 (r) ∂KP ∂CP ID ( r ) u1 (r) ∂K ∂CP ID (r) u1 (r) ∂Td + + + + + + ∂CP ID (r)y 2 (r) ∂KP ∂CP ID (r)y 2 (r) ∂Ti ∂CP ID (r)y 2 (r) ∂Td ∂CP ID (r)u2 (r) ∂KP ∂CP ID (r)u2 (r) ∂Ti ∂CP ID (r)u2 (r) ∂Td Da questi sviluppi, appare evidente che bisogna calcolare anche le derivate parziali di CP ID (s) e di CP I (s) rispetto ai parametri del controllore. Costante derivativa : ∂CP ID (r) ∂Td = =0 ∂CP I (r) ∂Td = − TK2Ps KP s (1+Td /N s)2 Costante Integrativa : ∂CP ID ∂Td ∂CP I (r) ∂Td i Costante proporzionale : ∂CP ID (r) ∂KP 204 i = − TK2Ps = Ti s+1 Ti s + ∂CP I (r) ∂KP = Ti s+1 Ti s Td s T 1+ Nd s di C.Bettini, A.Guiggiani e T.Lorenzetti 22.5 Calcolo dell’Hessiana di J 22.5. Calcolo dell’Hessiana di J Ci sono varie possibilità di scelta per la matrice Ri presente nella formula iterativa per l’aggiornamento dei parametri. La scelta più semplice è la matrice di identità, che ci dà la direzione negativa del gradiente. In alternativa si può usare: ˆ Ri = T " est 0 # " ∂y ∂y (r) est (r) ∂r ∂r #T " + lest # " ∂u ∂u (r) est (r) ∂r ∂r #T (22.30) dt che ci dà una stima non polarizzata Gauss-Newton della matrice Hessiana di J. Sostituendo : Ri = ˆT 0 est ∂y (r) ∂KP ∂y (r) ∂Ti ∂y (r) ∂Td est ∂y (r) ∂KP ∂y (r) ∂Ti ∂y (r) ∂Td T + λest ∂u (r) ∂KP ∂u (r) ∂Ti ∂u (r) ∂Td est ∂u (r) ∂KP ∂u (r) ∂Ti ∂u (r) ∂Td T dt Da ora omettiamo est[] per alleggerire la notazione: Ri = ˆT 0 ∂y 2 ∂K2P ∂y ∂KP 2∂Ti ∂y ∂KP ∂Td 2 ∂u ∂KP 2 λ ∂K∂u P ∂Ti λ ∂K∂u2 P ∂Td +λ + + 2 ∂y 2 + λ ∂K∂u ∂KP ∂Ti ∂T 2 P 2i ∂y ∂u + λ ∂KP2 T ∂Ti2 ∂y + λ ∂T∂u ∂Ti ∂Td i ∂Td ∂y 2 ∂KP ∂Td ∂y 2 ∂Ti ∂Td ∂y 2 ∂Td +λ + + ∂u2 ∂KP2∂Td λ ∂T∂u i ∂T 2d ∂u λ ∂T d Come si può notare è una matrice simmetrica quindi le componenti della matrice da calcolare si riducono. Passando al dominio di s si ottiene: 1 Ri = s ∂y 2 ∂K2P ∂y ∂KP 2∂Ti ∂y ∂KP ∂Td 2 ∂u ∂KP 2 λ ∂K∂u P ∂Ti λ ∂K∂u2 P ∂Td +λ + + 2 ∂y 2 + λ ∂K∂u ∂KP ∂Ti ∂Ti P 2 2 ∂y ∂u + λ ∂KP2 T ∂Ti2 ∂y + λ ∂T∂u ∂Ti ∂Td i ∂Td ∂y 2 ∂KP ∂Td ∂y 2 ∂Ti ∂Td ∂y 2 ∂Td di C.Bettini, A.Guiggiani e T.Lorenzetti +λ + + ∂u2 ∂KP2∂Td λ ∂T∂u i ∂T 2d ∂u λ ∂Td 205 dt 22.5 Calcolo dell’Hessiana di J (22.31) 206 di C.Bettini, A.Guiggiani e T.Lorenzetti 23. Costruzione dello schema Simulink Il processo in esame è controllato ed osservato tramite una scheda di acquisizione di tipo Real-Time che lavora su dati campionati. Quindi si è reso necessario discretizzare il metodo precedentemente trattato in tempo continuo. Per passare dal dominio di Laplace alla trasformata Zeta si può fare riferimento al seguente ragionamento: data una funzione di trasferimento F dT = Ti s + 1 y(s) = u(s) Ti s si può ricavare: y (Ti s) = u (Ti s + 1) e, quindi: y =u+u 1 sTi che riportato nel dominio Zeta equivale a 207 23.1 Schema Simulink y(z) = u + u 1 Ts Ti z − 1 23.1. Schema Simulink In Figura 23.1 si riporta lo schema Simulink del sistema di controllo dove possiamo notare i sottoblocchi temporizzati: • dal sottosistema “Primo_Esperimento” si ricavano i segnali di errore e uscita; • nel sottosistema “Secondo_Esperimento” si esegue il secondo esperimento e si determinano i parametri dell’Hessiana e il gradiente del costo; • fra i due blocchi sono presenti dei blocchi “delay” al fine di eseguire in sequenza temporale i due esperimenti riportando come riferimento del secondo esperimento l’errore del primo; (ri )) • il sottosistema “Algoritmo_di_aggiornamento” esegue il calcolo (gi Ri−1 ∂J ∂r necessario per l’aggiornamento dei parametri (vedi ??eq:algoritmo a pagina 196); • il sottosistema “Aggiornamento_parametri” aggiorna i parametri necessari all’iterazione e li mantiene fino all’iterazione successiva. Riportiamo di seguito in dettaglio i sottoblocchi appena descritti con le relative temporizzazioni. Il periodo di esecuzione di ogni iterezione è stato scelto pari a 12 secondi. 23.2. Primo Esperimento Il sottosistema “Primo_Esperimento” (Figura 23.2) è composto dal controllore PID collegato direttamente al sistema reale ed attivo nei primi 8s. 208 di C.Bettini, A.Guiggiani e T.Lorenzetti 23.3 Secondo Esperimento Figura 23.1.: Diagramma a blocchi completo del controllore Figura 23.2.: Schema del primo esperimento 23.3. Secondo Esperimento Il sottosistema “Secondo_Esperimento” (Figura 23.3) comprende l’esecuzione della seconda prova (con il controllore PID) ed il relativo calcolo dei parametri necessari all’algoritmo di aggiornamento. Il blocco è attivo tra gli 8s ed i 12s. Ricordiamo che in ingresso al secondo esperimento si trova l’errore (e = y−rif erimento) realtivo al primo esperimento e ritardato temporalmente di 8 secondi. di C.Bettini, A.Guiggiani e T.Lorenzetti 209 23.4 Algoritmo di aggiornamento Figura 23.3.: Schema del blocco Secondo Esperimento comprensivo del calcolo dei parametri 23.4. Algoritmo di aggiornamento Il blocco “Algoritmo_di_aggiornamento” è composto da un primo sottosistema che calcola la matrice inversa e da un parallelo di tre sottosistemi che calcolano gli incrementi (4K, 4Ti , 4Td ) dei parametri del PID (Figura 23.4). Figura 23.4.: Schema simulink del blocco “Algoritmo di Aggiormento” per il calcolo degli incrementi 23.5. Aggiornamento parametri Il blocco “Aggiornamento_parametri”, infine, è strutturato in modo tale da aggiornare alla fine di ogni periodo (più precisamente un campione prima della fine dell’iterazione) i parametri sfruttando un trigger e mantenendoli per tutto il periodo successivo (Figura 23.5). Figura 23.5.: Schema Simulink per l’aggiornamento parametri 210 di C.Bettini, A.Guiggiani e T.Lorenzetti 23.6 Segnale di controllo 23.6. Segnale di controllo Per quanto riguarda il segnale di controllo, facendo riferimento al diagramma a blocchi illustrato in Figura 23.1, si nota che il primo esperimento genera un segnale di controllo per i primi 8 secondi del ciclo mentre il secondo ne genera uno per i restanti 4 secondi. Oltre a ciò, l’errore in uscita dal primo esperimento è stato temporalmente filtrato (sfruttando le proprietà dei sottosistemi “enable”) in modo da avere in ingresso al secondo esperimento solo la parte riguardante l’errore relativo ai primi 4 secondi per evitare di introdurre segnali di errore non necessari per la corretta esecuzione dell’algoritmo. Inoltre il segnale di errore è ridotto di un fattore 1/2 per ovviare a problemi di esecuzione dell’esperimento. di C.Bettini, A.Guiggiani e T.Lorenzetti 211 24. Prove con riferimento a onda quadra In questo capitolo andremo ad analizzarere il metodo iterativo IFT. Per analizzare l’algoritmo abbiamo deciso di applicare diversi segnali di riferimento al sistema variando i parametri inizializzati. Inizialmente il sistema dovrà inseguire un riferimento ad onda quadra; in seguito saranno implementati anche riferimenti a dente di sega e sinusoidali. E’ stata infine analizzata la robustezza del metodo IFT quando il sistema è sottoposto a disturbi esterni. Il riferimento è stato scelto come un’onda quadra avente periodo 12 e un duty-cicle di 1/3. Per quanto riguarda l’errore, il diagramma Simulink è stato temporizzato in modo che esso fosse introdotto dopo 8 secondi (sempre con periodo pari a 12 secondi). Come primo esperimento è stato scelto un insieme di parametri che fosse sub-ottimo per quanto riguarda le prestazioni del sistema. Dopo ciò, è stata testata l’efficienza dell’algoritmo allontanandoci dalla zona di ottimalità dei parametri. E’ possibile osservare in Figura 24.1 l’andamento dell’uscita e dei parametri, inizializzati K = 1, T i = 1, T d = 0.04. Come si può notare, ad ogni passo il tempo di assestamento si riduce a scapito però di un aumento della sovraelongazione. In questo esempio il passo di iterazione è stato scelto pari a 0.04. Andiamo adesso a vedere cosa accade se viene incrementato del 50% il valore della costante derivativa (T d = 0.06). 213 Prove con riferimento a onda quadra Figura 24.1.: Uscita e parametri nel caso in cui K=1, Td=0.04, Ti=1, passo=0.04. Come si può notare in Figura 24.2 e Figura 24.3 l’uscita del sistema ed il segnale di controllo sono molto oscillanti anche se si nota comunque un miglioramento ad ogni iterazione. Figura 24.2.: Uscita nel caso in cui K=1, Td=0.06, Ti=1, passo=0.04. La simulazione appena descritta potrebbe far credere che l’algoritmo non lavori in modo corretto; incrementando però il passo, più precisamente raddoppiandolo, si nota un netto miglioramento rappresentato in Figura 24.4. Come si osserva in Figura 24.4, un passo più grande fa sì che si giunga abbastanza 214 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove con riferimento a onda quadra Figura 24.3.: Controllo nel caso in cui K=1, Td=0.06, Ti=1, passo=0.04. rapidamente in una zona in cui i parametri sono ottimi per poi allontanarci da essa a causa della struttura dell’algoritmo che non riesce a convergere sul minimo locale essendo il passo troppo elevato (Figura 24.5). Andiamo ora ad analizzare in Figura 24.6 il comportamento dell’uscita se i parametri sono inizializzati a K = 1, T i = 0.3, T d = 0.04 e con un passo = 0.06. Anche in questo caso, ad ogni iterazione si assiste ad una diminuzione dell’errore (sia in ampiezza che in durata) che viene riportato sull’uscita dopo i primi 8 secondi e, successivamente, con periodo 12 secondi. di C.Bettini, A.Guiggiani e T.Lorenzetti 215 Prove con riferimento a onda quadra Figura 24.4.: Uscita e segnale di controllo nel caso in cui K=1, Td=0.06, Ti=1, passo=0.08. Figura 24.5.: Parametri relativi al caso precedente 216 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove con riferimento a onda quadra Figura 24.6.: Uscita nel caso in cui K=1, Ti=0.3, Td=0.04, passo=0.06. di C.Bettini, A.Guiggiani e T.Lorenzetti 217 25. Prove con riferimento a dente di sega Analizziamo adesso il metodo IFT con un riferimento a dente di sega che varia da −1V a −3V con periodo 6 secondi. Dato che - per nostra scelta - la temporizzazione dello schema Simulink è invariata, ne consegue che l’errore verrà applicato una volta ogni due segnali. Osserviamo in Figura 25.1 il segnale d’uscita e l’andamento dei parametri (inizializzati K = 1, T d = 0.04, T i = 0.4 con passo = 0.03) mentre nella Figura 25.2 è raffigurato il segnale di controllo. Anche in questo caso si assiste ad un miglioramento e ad un affinamento del controllo che si manifesta nella riduzione dell’errore sia in ampiezza che in durata. Proviamo infine a rilassare l’esperimento precedente andando a variare Ti , più precisamente inizializzando K = 1, T d = 0.04, T i = 0.1 e ponendo il passo = 0.03. Oltre a ciò è stato aumentato il periodo del riferimento a dente di sega portandolo a 12 secondi in modo da avere una iterazione ogni ciclo. Si noti in particolare in Figura 25.4 quanto il segnale di controllo sia più moderato rispetto a quello mostrato in Figura 25.2. 219 Prove con riferimento a dente di sega Figura 25.1.: Uscita e parametri del sistema con riferimento a dente di sega e K=1, Td=0.04, Ti=0.4, passo=0.03. Figura 25.2.: Controllo ed il relativo ingrandimento 220 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove con riferimento a dente di sega Figura 25.3.: Uscita con riferimento a dente di sega e K=1, Td=0.04, Ti=0.1, passo=0.03. Figura 25.4.: Segnale di controllo relativo all’esperimento precedente. di C.Bettini, A.Guiggiani e T.Lorenzetti 221 26. Prove con riferimento sinusoidale Analizziamo infine l’efficacia del metodo finora descritto con un riferimento sinusoidale. Affinchè l’errore sia indotto in una parte del riferimento tale da eseguire correttamente gli esperimenti è stata scelta una sinusoide con periodo pari a 12 secondi ed ampiezza variabile tra −1V e −3V . Come vederemo all’interno di questo paragrafo, nel caso in cui il riferimento da inseguire è una sinusoide l’ottimizzazione è più efficiente in quanto questo tipo di segnali sono più compatibili con i regolatori PID. Osserviamo in Figura 26.1 l’andamento dell’uscita inizializzando i parametri a K = 1, T d = 0.04, T i = 0.17 e passo = 0.08. Analizziamo ora il caso in cui la costante integrativa viene aumentata portandola a T i = 0.8 (Figura 26.3). Per motivi di convergenza è stato dimezzato il passo di iterazione rispetto al caso precedente (passo = 0.04). Quest’ultimo caso è la prova che, pur partendo da parametri iniziali non buoni (quasi al limite dell’instabilità), il metodo iterativo IFT riesce comunque ad ottimizzare in un numero limitato di passi il controllo del sistema. 223 Prove con riferimento sinusoidale Figura 26.1.: Caso con ingresso sinusoidale e parametri K=1, Td=0.04, Ti=0.17, passo=0.08. Figura 26.2.: Parametri relativi al caso in cui l’ingresso è sinusoidale e vengono inizializzati K=1, Td=0.04, Ti=0.8, passo=0.04. 224 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove con riferimento sinusoidale Figura 26.3.: Caso con ingresso sinusoidale e parametri K=1, Td=0.04, Ti=0.8, passo=0.04. Figura 26.4.: Segnale di controllo e parametri relativi all’esperimento precedente. di C.Bettini, A.Guiggiani e T.Lorenzetti 225 27. Robustezza del metodo In questo paragrafo andremo ad analizzare la robustezza del metodo introducendo un disturbo esterno nel sistema durante l’esecuzione dell’algoritmo iterativo. Effettuiamo un primo test durante il quale il disturbo è indotto spostando la pallina con una penna. Si noti che, nel caso il riferimento sia di tipo sinusoidale come in Figura 27.1 e Figura 27.2, il sistema torna a seguire il riferimento non appena il disturbo cessa. Figura 27.1.: Caso in cui viene introdotto un disturbo esterno con riferimento sinusoidale. Vediamo ora in Figura 27.3 il comportamento del sistema nel caso in cui il disturbo sia indotto quando il riferimento è a dente di sega. 227 Robustezza del metodo Figura 27.2.: Caso in cui viene introdotto un disturbo esterno marcato con riferimento sinusoidale. Anche in questo caso il metodo iterativo si dimostra robusto a disturbi esterni. Analizziamo infine il caso in cui il segnale da inseguire è un’onda quadra: anche stavolta, come si nota in Figura 27.4, il disturbo viene reiettato. Andiamo adesso ad analizzare il comportamento del sistema nel caso in cui la pallina venga tolta dalla sua sede e quindi il disturbo si presenterà come molto accentuato. Come si può notare in Figura 27.5, il disturbo non intacca la robustezza del sistema sia nel caso il riferimento sia un’onda quadra, che nel caso in cui il riferimento sia una sinusoide. In entrambi i casi, infatti, una volta rimessa la pallina nella sede il sistema torna stabile e continua ad ottimizzare i parametri. 228 di C.Bettini, A.Guiggiani e T.Lorenzetti Robustezza del metodo Figura 27.3.: Caso in cui viene introdotto un errore con riferimento a dente di sega. Figura 27.4.: Caso in cui viene introdotto un disturbo esterno con riferimento a onda quadra. di C.Bettini, A.Guiggiani e T.Lorenzetti 229 Robustezza del metodo Figura 27.5.: Caso in cui viene introdotto un disturbo esterno con riferimento a onda quadra. 230 di C.Bettini, A.Guiggiani e T.Lorenzetti 28. Conclusioni In questo capitolo sono illustrate le conclusioni derivanti dal lavoro svolto in laboratorio. Di seguito sono mostrati i segnali di controllo, riferimento ed uscita del sistema una volta implementati i PID aventi i parametri ottimi ricavati grazie al metodo IFT. In ognuno dei seguenti esperimenti è stata implementata una rampa necessaria a portare la pallina dalla posizione iniziale (2.3V ) all’inseguimento del segnale di riferimento. Nel caso in Figura 28.1 è raffigurato l’inseguimento di un’onda quadra implementando il PID con i parametri ottimi K = 1.1298, T i = 0.4328, T d = 0.0409. In Figura 28.2, invece, è stato implementato un PID caratterizzato dai parametri ottimi K = 1.0458, T i = 0.3545, T d = 0.0398 nel caso di riferimento a dente di sega. Infine in Figura 28.3 si può osservare l’inseguimento di una sinusoide implementato un PID avente i parametri ottimi K = 1.5171, T i = 0.1365, T d = 0.0369. In conclusione possiamo affermare che il metodo IFT si è dimostrato efficace nella taratura di un PID ottimo pur non avendo una modellazione esatta del sistema. La potenzialità di questo metodo, infatti, non essendo model based, è il riuscire ad ottimizzare un controllore che di per sè stabilizza il sistema senza aver conoscenza del sistema stesso. 231 Conclusioni Figura 28.1.: Segnale di controllo e di uscita implementando il PID ottimo per riferimento a onda quadra. Figura 28.2.: Segnale di controllo e di uscita implementando il PID ottimo per riferimento a dente di sega. 232 di C.Bettini, A.Guiggiani e T.Lorenzetti Conclusioni Figura 28.3.: Segnale di controllo e di uscita implementando il PID ottimo per riferimento sinusoidale. di C.Bettini, A.Guiggiani e T.Lorenzetti 233 Parte VI. Laboratorio di automatica II 235 29. Controllo di un processo termico PT326 tramite PID con logica Fuzzy tarato attraverso algoritmi genetici 29.1. Introduzione All’interno di questa relazione è esaminato l’inseguimento del set point da parte di un processo termico (PT326) controllato con un PID la cui uscita è pesata con regole Fuzzy. Questa scelta è stata presa in quanto i PID che implementano la logica fuzzy riescono a controllare una più ampia gamma di sistemi rispetto ai classici PID industriali, ad esempio sistemi non lineari, di ordine elevato e con ritardo consistente. La potenza del fuzzy consta nel fatto che, come vedremo, non applicando un controllo lineare si riesce ad ottenere migliori prestazioni rispetto al caso classico. Il controllore esaminato ha una struttura PI+D in cui l’azione derivativa agisce sull’uscita e non sull’errore. Si ha così un doppio segnale di controllo (uPI e uD) il quale, prima di essere posto in ingresso al sistema, viene pesato da due blocchi fuzzy. Lo scopo principale del presente elaborato non è tanto il controllo di tale processo il quale, come vedremo, presenta problematiche di ritardo risolvibili solamente tramite l’uso di predittori, quanto lo studio di algoritmi genetici attraverso cui ottimizzare i sette parametri del sistema. Come vedremo, infatti, il controllore PI+D pesato fuzzy è abbastanza complesso e caratterizzato da sette parametri difficili da accordare tra loro. Per far ciò si è ricorso ad algoritmi genetici, metodi di ottimizzazione che si basano sul concetto di evoluzione e mutazione. Questi non sono algoritmi che 237 29.2 Descrizione del Processo ricercano un ottimo basandosi sulla minimizzazione di una funzione di costo ma creando una popolazione (scelta almeno all’inizio in modo casuale) che si evolve in modo tale da essere composta sempre da individui migliori (ovvero da individui che, nel nostro caso, migliorano le prestazioni del sistema) senza convergere (grazie alle operazioni di crossover e mutazione). È questa la differenza fondamentale tra algoritmi classici, che convergono ad un certo punto di minimo (locale o globale), e genetici, che si evolvono verso soluzioni sempre migliori esplorando in modo più completo lo spazio delle soluzioni. Una volta identificato l’impianto e tarato il controllo con gli algoritmi genetici in Matlab, andremo ad implementare fisicamente i parametri ottimi sul processo (controllabile tramite un VI Labview) notando infine come le risposte dei due ambienti di simulazione siano molto simili tra loro. 29.2. Descrizione del Processo Il processo termico PT326 modella situazioni in cui si ha la necessità di controllare la temperatura di un condotto (o più in generale di un ambiente) in presenza di disturbi esterni e di una funzione di trasferimento con ritardo. Il processo (illustrato in figura) consiste in un ventilatore centrifugo che preleva aria dall’esterno la quale viene riscaldata da una resistenza e poi incanalata attraverso un condotto al termine del quale è presente una termocoppia (Bead Thermistor). La velocità del ventilatore è costante mentre è possibile variare l’apertura del parzializzatore (Throttle) da cui si preleva l’aria. Ci troveremo così di fronte a diverse funzioni di trasferimento caratterizzate da risposte differenti. Dal punto di vista di ingresso ed uscita il sistema è rappresentabile come: 238 di C.Bettini, A.Guiggiani e T.Lorenzetti 29.2 Descrizione del Processo Figura 29.1.: Processo termico PT326 Figura 29.2.: Schematizzazione ingresso-uscita del processo dove n(t) caratterizza il disturbo, u(t) l’ingresso ed y(t) l’uscita misurata dalla termocoppia. A seguito di prove empiriche si dimostra che la funzione di trasferimento del processo termico tra ingresso e uscita è parametrizzabile nel dominio di Laplace di C.Bettini, A.Guiggiani e T.Lorenzetti 239 29.3 Identificazione del Processo come una funzione di trasferimento del primo ordine con ritardo: P (s) = K e−τ s s−p Dove il guadagno K, il polo p e il ritardo t sono variabili a seconda dell’apertura del parzializzatore e della temperatura esterna. Per quanto riguarda le specifiche di ingresso, il segnale di controllo è limitato ad essere compreso tra 0V e 10V . Nel nostro caso siamo andati ad identificare la funzione di trasferimento del processo con parzializzatore aperto rispettivamente a 30°, 60°, 90°, 100°, 120°, 160°. 29.3. Identificazione del Processo Per identificare il processo è stato anzitutto necessario raccogliere i dati riguardanti la sua risposta. Nel caso in esame è stato creato un VI in Labview in cui in ingresso al processo è applicata un’onda quadra. I dati acquisiti sono stati poi salvati attraverso un secondo VI in un file di testo e rielaborati tramite Matlab sfruttando il metodo del simplesso flessibile per poter così identificare al meglio il processo. Con la tecnica del simplesso flessibile si cerca di minimizzare un funzionale d’errore: E(P ) = Nexp 1 X kYexp − Ym k2 Nexp i=1 Ricercando una stima: P̄ = arg min{E(P )} . P̄ ∈ RNp P Il simplesso flessibile è un metodo di ricerca diretta basato sulla conoscenza per punti del funzionale d’errore E(P). Questo viene valutato negli n+1 vertici del simplesso (dove n è il numero di parametri) e, ad ogni passo, viene scartato l’elemento peggiore (ovvero quello con E(P) massimo). È da notare che tutti gli elementi sono pesati con probabilità uniforme. L’elemento selezionato viene quindi sostituito da uno migliore 240 di C.Bettini, A.Guiggiani e T.Lorenzetti 29.3 Identificazione del Processo sfruttando l’algoritmo di Fibonacci in modo da localizzare i minimi locali per il funzionale d’errore. Il metodo si considera concluso quando tutti i valori di E(P) sono simili. Usando i dati ricavati dal processo fisico (avente in ingresso un’onda quadra di ampiezza 3) ed implementando il tutto nell’algoritmo si ottengono i risultati illustrati nelle figure seguenti e ricapitolati nella tabella successiva. Durante la fase sperimentale è stato notato che i dati raccolti sono influenzati anche dalla temperatura esterna. Come vedremo è risultato difficile ottenere un modello preciso del processo poiché esso è soggetto a numerosi disturbi esterni. Figura 29.3.: Modello con parzializzatore a 30° Figura 29.4.: Modello con parzializzatore a 60° di C.Bettini, A.Guiggiani e T.Lorenzetti 241 29.3 Identificazione del Processo Figura 29.5.: Modello con parzializzatore a 90° Figura 29.6.: Modello con parzializzatore a 100° Dato il modello in forma: P (s) = K e−sτ + of f set s+p A questo punto è stato necessario discretizzare la funzione di trasferimento del processo. Per fare ciò è stato utilizzato il comando Matlab “c2d” che fornisce direttamente il modello in tempo discreto (dato un certo tempo di campionamento che nel nostro caso è posto pari a 0.1 secondi). Ne derivano le funzioni di trasferimento, rispettivamente a 30°, 60°, 90°, 100°, 120° e 160°: G30° (z) = z −4 242 0.1732 (z − 0.7655) di C.Bettini, A.Guiggiani e T.Lorenzetti 29.3 Identificazione del Processo Figura 29.7.: Modello con parzializzatore a 120° Figura 29.8.: Modello con parzializzatore a 160° G360° (z) = z −3 G90° (z) = z −3 0.138 (z − 0.7338) 0.1129 (z − 0.7291) G100° (z) = z −3 0.1105 (z − 0.7183) G120° (z) = z −2 0.09483 (z − 0.7305) G160° (z) = z −2 0.09368 (z − 0.698) Si noti che il tempo di campionamento scelto è di 0.1 secondi, il quale è risultato un di C.Bettini, A.Guiggiani e T.Lorenzetti 243 29.3 Identificazione del Processo Figura 29.9.: Dati relativi ai modelli considerate le diverse parzializzazioni buon compromesso tra dinamica del processo e disturbi. 244 di C.Bettini, A.Guiggiani e T.Lorenzetti 30. PI+D con logica fuzzy E’ noto che i PID convenzionali non lavorano bene con sistemi non lineari, di ordine elevato o con ritardi incisivi. Per far sì che queste caratteristiche siano migliorate, è possibile usare un controllore PID adattativo impiegando una logica fuzzy. Questo controllore fuzzy ha delle caratteristiche interessanti in quanto: 1. Ha la stessa struttura lineare di un PID convenzionale ma con guadagni adattativi; i guadagni proporzionale, integrale e derivativo sono infatti funzioni non lineari del segnale di ingresso; 2. Il controllore è costruito sulla base di un PID discreto da cui vengono originate le leggi di controllo fuzzy; 3. Le funzioni-membro sono tutte funzioni triangolari e le regole sono solo quattro (tutte IF-THEN); ne deriva una legge di controllo semplice e comunque performante. Per quanto riguarda la stabilità di questi controllori PID, essa non è trattata all’interno di questa relazione. È però possibile fare riferimento alle letture in bibliografia per maggiori informazioni. Nel caso in esame il controllore scelto è del tipo PI+D, un controllore in cui l’azione derivativa viene applicata sull’uscita del processo e Figura 30.1.: Struttura di un controllore PI+D 245 30.1 Discretizzazione del PID 30.1. Discretizzazione del PID Andiamo ora a trattare la discretizzazione prima del controllore PI e poi del controllore D. 30.2. Controllore PI: L’uscita del controllore PI è esprimibile nel dominio di Laplace come: uP I (s) = Kpc Kc + i e(s) s Dove Kpc e Kic sono i guadagni, proporzionale ed integrativo, mentre e(s) è il segnale di errore. Questa equazione può essere riportata in tempo discreto applicando la trasformazione bilineare: 2 s= Ts z+1 z−1 in cui Ts è il tempo di campionamento. Ne deriva: uP I (z) = Kpc K c Ts Kic Ts − i + e(z) 2 1 − z −1 Imponendo: Kp = Kpc Kic Ts − 2 Ki = Kic Ts e passando dal dominio Zeta ai dati campionati, si ottiene: uP I (nTs ) − uP I (nTs − Ts ) = Kp [e(nTs ) − e(nTs − Ts )] + Ki Ts e(nTs ) 246 di C.Bettini, A.Guiggiani e T.Lorenzetti 30.3 Controllore D: da cui: uP I (nTs ) = uP I (nTs − Ts ) + Ts ∆uP I (nTs ) con: ∆uP I (nTs ) = Kp ev (nTs ) + Ki ep (nTs ) ev (nTs ) = e(nTs ) − e(nTs − Ts ) Ts ep (nTs ) = e(nTs ) Dove ∆uP I (nTs ) è il controllo incrementale del controllore PI, ep (nTs ) è il segnale di errore ed ev (nTs ) caratterizza la variazione dello stesso. A questo punto si può applicare una azione di controllo fuzzy al posto del termine . Si ottiene quindi: uP I (nTs ) = uP I (nTs − Ts ) + KuP I ∆uP I (nTs ) 30.3. Controllore D: Per quanto riguarda l’azione derivativa, il segnale di controllo è dato da: uD (s) = sKdc y(s) con Kdc pari al guadagno derivativo e y(s) al segnale di uscita. Usando nuovamente la trasformazione bilineare si ottiene: uD (z) = Kdc 2 1 − z −1 y(z) Ts 1 + z −1 da cui: uD (nTs ) + uD (nTs − Ts ) = 2Kdc [y(nTs ) − y(nTs − Ts )] Ts di C.Bettini, A.Guiggiani e T.Lorenzetti 247 30.3 Controllore D: Passando ai dati campionati: uD (nTs ) = −uD (nTs − Ts ) + KuD ∆uD (nTs ) Dove KuD è un guadagno costante del controllore che dovrà essere determinato e: ∆uD (nTs ) = Kd ∆y(nTs ) con: ∆y(nTs ) = y(nTS ) − y(nTs − Ts ) Ts pari alla variazione dell’uscita e: Kd = 2Kdc Ts Per ottenere migliori prestazioni del controllore derivativo è possibile modificare sommandoci il segnale yd (nTs ): yd (nTs ) = y(nTs ) − r(nTs ) = −e(nTs ) così da ottenere: ∆uD (nTs ) = Kd ∆y(nTs ) + Kd y(nTs ) Si noti che il coefficiente K serve per poter pesare anche yd (nTs ) in ingresso al regolatore Fuzzy. 248 di C.Bettini, A.Guiggiani e T.Lorenzetti 30.4 Controllore PI+D: 30.4. Controllore PI+D: Sommando il controllo PI al fuzzy D ne deriva il controllore PI+D che ha seguente legge di controllo: uP ID (nTs ) = uP I (nTs ) + KuP I ∆uP I (nTs ) + uD (nTs − Ts ) − KuD ∆uD (nTs ) di C.Bettini, A.Guiggiani e T.Lorenzetti 249 30.4 Controllore PI+D: Figura 30.2.: Schema del controllore PI+D con logica fuzzy 250 di C.Bettini, A.Guiggiani e T.Lorenzetti 31. Controllo Fuzzy Andiamo ora ad analizzare i passi per costruire il modello Fuzzy. Figura 31.1.: Schema generico del sistema inferenziale fuzzy 31.1. Fuzzificazione: Le componenti PI e D del controllore PI+D vengono fuzzificate in modo individuale. Per quanto riguarda il controllore fuzzy PI, esso è caratterizzato da due ingressi: aaa aaa Allo stesso modo il controllore fuzzy D ha in ingresso: bbb 251 31.1 Fuzzificazione: bbb Per quanto riguarda le uscite, il controllore PI genera un segnale mentre il controllore D genera un segnale . Come si può notare dalle regole Fuzzy osservabili in Figura 31.2 e Figura 31.3, viene utilizzata una sola costante (L) per le funzioni d’appartenenza. Gli ingressi sono pesati con le costanti Kp , Ki , Kd e K mentre le uscite con le costanti Ku P I e Ku D. Figura 31.2.: Funzioni membro del componente PI Figura 31.3.: Funzioni membro del componente D 252 di C.Bettini, A.Guiggiani e T.Lorenzetti 31.2 Regole d’inferenza: 31.2. Regole d’inferenza: Le regole d’inferenza fuzzy per il controllore PI sono quindi definite come: Allo stesso modo per il controllore D si ha che: IF ep = negativo AN D ev = negativo T HEN P Ioutput = negativo IF ep = negativo AN D ev = positivo T HEN P Ioutput = zero IF ep = positivo AN D ev = negativo T HEN P Ioutput = zero IF ep = positivo AN D ev = positivo T HEN P Ioutput = positivo Allo stesso modo per il controllore D si ha che: IF ỹd = negativo AN D ∆ỹ = negativo T HEN P Ioutput = zero IF ỹd = negativo AN D ∆ỹ = positivo T HEN P Ioutput = positivo IF ỹd = positivo AN D ∆ỹ = negativo T HEN P Ioutput = negativo IF ỹd = positivo AN D ∆ỹ = positivo T HEN P Ioutput = zero Per capire meglio le regole appena descritte vediamo come viene interpretata la prima regola del controllore PI. Se l’errore e la sua variazione sono negativi allora l’uscita (y) si trova al di sopra del setpoint e sta crescendo (perché l’errore a derivata negativa implica l’uscita con derivata positiva dato che e = ysp −y ed ysp è costante); quindi l’azione ∆uP I (nTs ), vista nelle formule precedenti, che riguarda il controllo PI contenente l’integratore viene posta negativa mentre la componente derivativa ∆uD (nTs ) è zero. 31.3. Defuzzificazione Il concetto impiegato per defuzzificare gli incrementi dei controllori PI e D è descritto dalla formula: di C.Bettini, A.Guiggiani e T.Lorenzetti 253 31.3 Defuzzificazione ∆u(t) = P (valori d0 ingresso × uscita corrispondente al valore d0 ingresso) P valori d0 ingresso Mettendo insieme le funzioni di appartenenza di ẽp ed ẽv rispettivamente sull’asse delle ascisse e sull’asse delle ordinate si ottiene una funzione in tre dimensioni dove le proiezioni della terza dimensione sul piano x, y creano 20 regioni. Prendendo ad esempio la regione IC1 (in Figura 31.4) possiamo notare i domini ep ∈ [0, L] ed ev ∈ [0, L]. L’insieme tra le regole d’inferenza del controllore Fuzzy con le funzioni d’appartenenza e le regioni serve a valutare le appropriate leggi di defuzzificazione. Quindi considerando ed nella regione IC1 notiamo che: ẽv negativo > 0.5 > ẽp negativo Facendo riferimento alla regola 1 del fuzzy PI (considerando che l’operatore AND prende il minimo) si ottiene: regola1 = se l0 ingresso selezionato (ẽp ) e0 negativo allora il corrispondente valore d0 uscita e0 negativo Seguendo lo stesso ragionamento per le altre tre regole si ottiene: regola2 = regola3 = 254 se l0 ingresso selezionato (ẽp ) e0 positivo allora il corrispondente valore d0 uscita e0 zero se l0 ingresso selezionato (ẽv ) e0 negativo allora il corrispondente valore d0 uscita e0 zero di C.Bettini, A.Guiggiani e T.Lorenzetti 31.3 Defuzzificazione se l0 ingresso selezionato (ẽv ) e0 positivo allora regola4 = il corrispondente valore d0 uscita e0 positivo Si può verificare che queste affermazioni sono valide anche per la regione IC2; applicando la formula di defuzzificazione si ottiene: ∆u(t) = ẽp negativo × outn + ẽp negativo × outz + ẽv negativo × outz + ẽv positivo × outp ẽp negativo + ẽp negativo + ẽv negativo + ẽv positivo In cui outn = −L, outz = 0 ed outp = L corrispondono ai relativi valori d’uscita. Figura 31.4.: Defuzzificazione del controllo PI (a) e del controllo D (b) Data la struttura geometrica delle funzioni membro (illustrata in Figura 31.2 e Figura 31.3) si ottiene: ẽp positivo = Ki ep (nTs ) + L 2L ẽp negativo = −Ki ep (nTs ) + L 2L ẽv positivo = Kp ev (nTs ) + L 2L di C.Bettini, A.Guiggiani e T.Lorenzetti 255 31.3 Defuzzificazione ẽv positivo = −Kp ev (nTs ) + L 2L Quindi, riferendosi alla formula precedente che definiva , si ottiene: ∆uP I (nTs ) = L [Ki ep (nTs ) + Kp ev (nTs )] 2 (2L − Ki ep (nTs )) Si può notare che ep (nTs ) > 0 è verificata per le regioni IC1 e IC2 mentre per le regioni IC5 e IC6 si ottiene: ∆uP I (nTs ) = L [Ki ep (nTs ) + Kp ev (nTs )] 2 (2L + Ki ep (nTs )) con ep (nTs ) < 0. Quindi è possibile adottare un’unica formula di defuzzificazione per queste quattro regioni (IC1,2,5,6): uP I (nTs ) = L [Ki ep (nTs ) + Kp ev (nTs )] 2 (2L + Ki |ep (nTs )|) Seguendo gli stessi ragionamenti per le altre regioni si ottengono le regole (illustrate in seguito) di defuzzificazione che generano le due uscite ∆uP I (nTs ) e ∆uD (nTs ). ∆uP I (nT ) = 256 L[Ki ep (nT )+Kp ev (nT )] in IC1, 2, 5, 6 2[2L−Ki ep (nT )] L[Ki ep (nT )+Kp ev (nT )] in IC3, 4, 7, 8 2[2L−Kp ev (nT )] 1 [Kp ev (nT ) + L] in IC9, 10 2 1 [Ki ep (nT ) + L] in IC11, 12 2 1 [Kp ev (nT ) − L] in IC13, 14 2 1 [Ki ep (nT ) − L] in IC15, 16 2 L in IC17 −L in IC19 di C.Bettini, A.Guiggiani e T.Lorenzetti 31.3 Defuzzificazione ∆uP I (nT ) = 0 in IC18, 20 ∆uD (nT ) = L[Kyd (nT )+Kd ∆y(nT )] in IC1, 2, 5, 6 2[2L−Kyd (nT )] L[Kyd (nT )+Kd ∆y(nT )] in IC3, 4, 7, 8 2[2L−Kd ∆y(nT )] 1 [−Kd ∆y(nT ) + L] in IC9, 10 2 1 [yd (nT ) + L] in IC11, 12 2 1 [−Kd ∆y(nT ) − L] in IC13, 14 2 1 [yd (nT ) − L] in IC15, 16 2 −L in IC18 L in IC20 ∆uD (nT ) = 0 in IC17, 19 di C.Bettini, A.Guiggiani e T.Lorenzetti 257 32. Algoritmi genetici 32.1. Introduzione agli algoritmi genetici Gli algoritmi genetici (GA) sono algoritmi di ricerca basati sui principi genetici della selezione naturale e dell’evoluzione: in essi si fondono i concetti di sopravvivenza degli elementi migliori, di riproduzione e di mutazione. In modo iterativo, partendo da una popolazione iniziale, un algoritmo genetico genera una nuova popolazione utilizzando gli elementi migliori della precedente ed applicandoci delle variazioni. Ad ogni elemento della popolazione viene poi associata una misura “prestazionale” relativa ad una funzione di fitness in modo tale che, alla fine di ogni ciclo, le stringhe con più alta affinità (fitness) si possono riprodurre. Gli algoritmi genetici per il loro modo di operare forniscono così un compromesso fra prestazioni ed esplorazione: ogni popolazione va a migliorare (sopprimendo gli esemplari più deboli) da un lato generando però in modo pseudocasuale (tramite la tecnica di crossover e la mutazione) dei figli che potrebbero portare (o meno) ad un miglioramento della specie. Le differenze fondamentali tra algoritmi genetici e metodi di ottimizzazione tradizionale sono riconducibili al fatto che gli algoritmi genetici: 1. lavorano con parametri codificati; 2. compiono la ricerca su una popolazione di punti e non considerando un dato alla volta; 3. utilizzano direttamente il valore della funzione oggetto (fitness) senza ricorrere ad informazioni aggiuntive come segni o derivate; 4. usano regole di transizione probabilistiche e non deterministiche. 259 32.2 Definizione del problema Oltre a ciò, negli algoritmi genetici non si va alla ricerca di un valore di minimo o di massimo come nel caso dell’ottimizzazione tradizionale: un GA non corre il rischio di fermarsi su un minimo locale ma, anzi, continua ad esplorare zone nuove e diverse (sebbene col rischio di andare incontro a soluzioni non ottimali). Quando si trattano gli algoritmi genetici è bene tenere a mente alcune definizioni: • Cromosoma: stringa, elemento della popolazione; • Gene: bit, elemento della stringa; • Locus: posizione della stringa; • Fitness: valore della funzione oggetto; • Genotipo: popolazione, insieme di tutti gli esemplari; • Fenotipo: parametro decodificato. 32.2. Definizione del problema Per poter applicare un algoritmo genetico occorre anzitutto codificare la popolazione; in genere la codifica adottata è la rappresentazione binaria, quindi ogni gene ha un valore 0 o 1. In seguito si dovrà definire un’opportuna funzione di fitness. 32.3. Operazioni genetiche Definite la popolazione, la codifica e la funzione di fitness, è possibile far evolvere la popolazione tramite gli operatori fondamentali degli algoritmi genetici: 1. Riproduzione 2. Crossover 3. Mutazione Durante la riproduzione ogni stringa viene duplicata in accordo con il suo valore della funzione oggetto (fitness) per generare la popolazione successiva: stringhe con valori alti di fitness hanno probabilità maggiore di riprodursi. 260 di C.Bettini, A.Guiggiani e T.Lorenzetti 32.3 Operazioni genetiche Il crossover, invece, procede in due passi: nel primo vengono accoppiati due cromosomi della popolazione, nel secondo si seleziona con probabilità uniforme un intero k compreso nell’intervallo [1, L] (dove L è la lunghezza della stringa) e si scambiano i geni tra i due genitori compresi nell’intervallo [k+1, L]. Ad esempio se si hanno due cromosomi di lunghezza L=5 con k=3, l’operazione di crossover avviene nel seguente modo: Prima del crossover: A = 100|10 B = 110|11 Dopo il crossover: A0 = 1 0 0 | 1 1 B0 = 1 1 0 | 1 0 Figura 32.1.: La mutazione, infine, consiste nel modificare casualmente il valore di uno o più bit della matrice contenente la popolazione. La presenza di questo operatore è giustificata dal fatto che l’alterazione casuale di un bit porta all’esplorazione di regioni completamente nuove dello spazio verso possibili soluzioni migliori. La frequenza con cui viene riscontrata la mutazione è stimata essere nell’ordine di una mutazione su mille bits. di C.Bettini, A.Guiggiani e T.Lorenzetti 261 32.4 Applicazioni per l’ottimizzazione del controllore 32.4. Applicazioni per l’ottimizzazione del controllore Dalla lettura del paragrafo precedente si evince quanto gli algoritmi genetici, oltre ad essere innovativi, riescano a lasciare un’ampia libertà di scelta per quanto riguarda funzioni di fitness e codifica dei dati. Nel nostro caso i parametri di controllo sono sette; questi devono innanzitutto essere codificati in una stringa binaria di n bits. Il cromosoma è così definito: I = (Kp Ki Kd L KuP I KuD K) La codifica è fatta facendo convertendo in binario tali parametri, specificando una precisione ed una risoluzione in bit uguale per tutti. La funzione di fitness è stata scelta in modo tale da minimizzare una serie di parametri. Tale funzione è una combinazione lineare dei parametri da minimizzare pesati opportunamente. I parametri scelti sono: • il valore di overshoot (os); • il tempo di salita (tr); • il tempo di assestamento (ts); • il valore di undershoot (us); • il valore dell’energia dell’errore (en). F = [os tr ts us en] L’energia dell’errore viene calcolata tramite la formula: ˆT en = [y(t) − ysp (t)]2 dt 0 262 di C.Bettini, A.Guiggiani e T.Lorenzetti 32.5 Implementazione Quindi la funzione di fitness diventa: f = w1 os + w2 tr + w3 ts + w4 en + w5 us Dove i pesi possono essere assegnati in modo casuale o preimpostato (per la trattazione dei pesi si rimanda al capitolo successivo). A questo punto si vanno a eseguire le simulazioni riportando in decimale ogni cromosoma e ricavando i dati per calcolare i parametri di fitness (opportunatamente pesati) a cui verranno associati. L’operazione di riproduzione avviene mediante la tecnica della roulette polarizzata. Con questa tecnica si associa ad ogni esemplare una certa probabilità a seconda del valore di fitness: in particolare nel nostro caso viene assegnata una probabilità maggiore per la riproduzione agli elementi che presentano un valore di fitness minore (e, quindi, sono caratterizzati da migliori prestazioni). In questo modo si crea una nuova popolazione riproducendo i cromosomi in base alla probabilità assegnata. L’operazione di crossover viene eseguita un numero di volte pari al 5% (arrotondato) del numero di esemplari appartenenti allla popolazione prendendo ciascuna volta una coppia casuale di cromosomi. La mutazione, infine, viene eseguita variando un bit dell’intera popolazione in maniera casuale con una frequenza di 1 su 1000 bits. 32.5. Implementazione In questo capitolo andremo a descrivere soluzioni ed accorgimenti messi in atto durante l’implementazione dell’algoritmo. È stata inizialmente generata una popolazione in modo casuale ed è stata fatta evolvere fino al raggiungimento di un esemplare che stabilizzasse il sistema. In un secondo momento da questo esemplare è stato codificato il primo cromosoma mentre gli altri sono stati codificati partendo dai parametri della soluzione stabilizzante, perturbandoli con un valore casuale al fine di generare una popolazione di C.Bettini, A.Guiggiani e T.Lorenzetti 263 32.5 Implementazione nell’intorno della zona di stabilità. Questo ci permette di arrivare ad una soluzione psedudo-ottima in un numero inferiore passi. Come detto nel capitolo precedente i pesi della funzione di fitness possono essere calcolati secondo un’assegnazione casuale o preimpostata: 1. Assegnazione casuale: I pesi dei parametri della fitness sono: osi, tri, tsi, eni, usi. Ad essi viene assegnato un valore casuale (attraverso la funzione rand) e vengono calcolati i coefficienti della funzione oggetto normalizzati come: w1 = osi ; osi+tri+tsi+en+usi w2 = tri ; osi+tri+tsi+en+usi w4 = en ; osi+tri+tsi+en+usi w5 = usi osi+tri+tsi+en+usi w3 = tsi osi+tri+tsi+en+usi 2. Assegnazione preimpostata: Nell’assegnazione preimpostata i pesi vengono scelti dall’utente, quindi: w1 = osi; w2 = tri; w3 = tsi; w4 = en; w5 = usi Tra le due opzioni la scelta è ricaduta sulla seconda in quanto in questo modo si ha un maggiore controllo sulle caratteristiche che si desidera minimizzare. Nel primo caso, invece, avendo i pesi casuali si vanno ad esplorare molteplici soluzioni non avendo una precisa direzione di ricerca. È da notare che per la ricerca della soluzione ottima non è stato impostato nessun criterio d’arresto ma semplicemente viene eseguito l’algoritmo in un ciclo for che si ferma dopo un definito numero di passi al fine di valutare la soluzione trovata ed, eventualmente, ri-inizializzarla per un nuovo ciclo. Per effettuare la simulazione del sistema di controllo si è costruito una schema Simulink, quindi è stato necessario caricare tutti i parametri riportandoli da binario 264 di C.Bettini, A.Guiggiani e T.Lorenzetti 32.5 Implementazione in decimale. Per poter modificare il punto di break-point L delle regole d’inferenza sono state usate le funzioni Matlab readfis e writefis. Un problema riscontrato è stato quello di ricavare i dati della risposta al gradino. Il sistema da controllare presenta un offset, quindi una simulazione Simulink che va a valutare le prestazioni sulla risposta ad un gradino unitario sarebbe errata. Per ovviare a ciò è stato implementato un generatore di onda quadra andando a valutare le prestazioni del sistema sul fronte di salita del secondo ciclo (così da evitare la parte caratterizzata dall’offset). Per poter calcolare l’energia dell’errore solo sul secondo fronte di salita abbiamo implementato un subsystem-enable con in ingresso il segnale d’errore che attiva l’uscita con un comando di enable opportunamente temporizzato. I segnali di tempo, uscita, riferimento ed energia dell’errore sono stati elaborati grazie ad un secondo m-file (RispostaGrad.m) da noi creato il quale calcola overshoot, tempo di salita, tempo di assestamento ed undershoot per ogni esemplare della popolazione. di C.Bettini, A.Guiggiani e T.Lorenzetti 265 32.5 Implementazione Figura 32.2.: 266 di C.Bettini, A.Guiggiani e T.Lorenzetti 33. Prove In questa parte finale della relazione andremo ad analizzare i risultati ottenuti applicando una serie di controllori ricavati attraverso gli algoritmi genetici appena descritti. E’ da notare che l’intero sistema è altamente variabile in funzione della temperatura esterna (e, inoltre, di un offset manuale). Si è quindi reso necessario modellare nuovamente il sistema ed effettuare le prove in un arco di tempo abbastanza ristretto poichè i modelli ricavati in precedenza non sarebbero stati esatti. Fortunatamente gli algoritmi appena descritti sono stati completamente automatizzati permettendoci di effettuare l’ottimizzazione in modo abbastanza rapido. Per quanto riguarda il modello, è stato preso in esame il sistema con parzializzatore a 100°. Ne è derivata una funzione di trasferimento nel dominio della frequenza: P (s) = 1.223 −0.25s e s + 3.169 La quale, riportata nel dominio Zeta, equivale a: P (z) = z2 0.1048 (z − 0.7284) Caratterizzata da un offset pari a 2.6091. A questo punto è stato fatto partire l’algoritmo inizializzando una popolazione di quelle ricavate nei capitoli precedenti. Ne è derivata una prima popolazione avente: 267 Prove Figura 33.1.: Modellazionedel processo con parzializzatore a 100° kp = 10 ki = 38, 7 kd = 0 kud = kupi = 0, 01 k = 1, 2 L = 58 La quale, implementata nell’algoritmo come popolazione iniziale, dopo cinque iterazioni da origine a: In questo esempio i pesi della fitness sono scelti in modo da minimizzare l’overshoot, risulta quindi che la risposta al gradino scelta come “migliore” sia la seguente: Il caso appena illustrato è caratterizzato da una tolleranza minima all’overshoot (pesando osi con un valore molto alto) a discapito di quella relativa al tempo di salita (che è pari a 0,9 secondi) avente peso minore. Se decidessimo quindi di cercare un controllo avente tempo di salita migliore a fronte di un piccolo overshoot, potremmo variare i parametri-peso (osi, tri, tsi, eni, usi) andando a diminuire quello relativo all’overshoot. 268 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove Figura 33.2.: Rappresentazione della popolazione (10 esemplari) dopo 5 cilci Ne deriva una diversa scelta della risposta migliore caratterizzata dai parametri: kp = 11 ki = 35 kd = 0 kud = kupi = 0, 01 k = 0, 8 L = 75 I quali danno origine ad una risposta al gradino con un tempo di salita migliore (circa 0,7 secondi) e un “leggerissimo” overshoot. A questo punto è stato possibile riportare tutto il controllo ottimizzato implementandolo in un VI Labview. Ne deriva lo schema in pagina seguente. A seguito di ulteriori iterazioni dell’algoritmo genetico sono stati ricavati dei parametri in cui nelle simulazioni Simulink era presente dell’overshoot: di C.Bettini, A.Guiggiani e T.Lorenzetti 269 Prove Figura 33.3.: migliore esemplare senza overshoot kp = 17, 9 ki = 71, 3 kd = 0, 01 kud = kupi = 0, 01 k=0 L = 157 La risposta del processo reale è molto simile a quella Simulink in Figura 33.7. caratterizzata da un overshoot del 16% e da un tempo di salita di 0,4 secondi metre il tempo di assestamento è pari a 1,4 secondi. Si noti inoltre come il processo sia di per sé abbastanza rumoroso e che il ritardo non può essere eliminato con questo tipo di controllo. A questo punto abbiamo cercato di ottimizzare un controllo che non presentasse nessun overshoot. A seguito dell’ottimizzazione genetica ne derivano i parametri: 270 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove Figura 33.4.: Migliore esemplare con leggero overshoot kp = 18.7 ki = 51, 3 kd = 0, 01 kud = kupi = 0, 01 k = 0, 01 L = 87 I quali danno una risposta in Simulink come quella in Figura 33.9 mentre la risposta del processo reale è illustrata in Figura 33.10. Si noti la differenza di risposta tra fronte di salita e fronte di discesa (in cui è presente un undershoot). Ciò è dovuto a una diversa risposta del sistema tra salita e discesa che potrebbe essere risolvibile tramite l’utilizzo di due modelli diversi (quindi controlli diversi). Il tempo di salita è pari a 0,8 secondi e il tempo di assestamento è minore rispetto al caso precedente (Figura 33.11). di C.Bettini, A.Guiggiani e T.Lorenzetti 271 Prove Se si implementa nuovamente l’algoritmo genetico si ottiene, dopo alcune iterazioni: kp = 11, 3 ki = 35, 9 kd = 0 kud = k = 0, 01 k = 0.5 L = 86 I quali in Simulink danno luogo alla risposta in Figura 33.12 che ha un corrispettivo nel processo reale con un andamento illustrato in Figura 33.13 ed in Figura 33.14. Notiamo in questa simulazione rispetto alla precedente l’assenza di undershoot e una migliore reiezione ai disturbi, il tempo di salita risulta essere pari a 0,8 secondi. Andiamo adesso a testare la robustezza di questo controllo: prendiamo i parametri ottimi dell’ultima simulazione (che sono anche quelli con il miglior compromesso overshoot-tempo di salita-tempo di assestamento) e andiamo a disturbare l’impianto aprendo e chiudendo il parzializzatore (prima aumentandolo fino a 160° e poi richiudendolo fino a 60°). Come possiamo notare in Figura 33.15 il sistema riesce comunque a controllare l’impianto. Notiamo inoltre in Figura 33.16 come al diminuire dell’apertura (quindi diminuendo l’afflusso d’aria) il segnale di controllo tende a decrescere. Proviamo adesso a chiudere ancora di più il parzializzatore (arrivando fino a 10°): ne deriva il risultato in Figura 33.17 con un segnale di controllo come quello in Figura 33.18. Come si può notare dalle figure, più il parzializzatore è chiuso più il controllo tende ad abbassarsi (data il minore afflusso di aria). Essendo il segnale di controllo costretto tra 0 e 10 V, si ha che nella parte finale del grafico il segnale di controllo (quando il parzializzatore è posto a 10°) porta in saturazione l’attuatore (Figura 33.19); non riuscendo ad inseguire il set point il segnale di controllo tende a divergere (poiché non è stata implementata alcuna tecnica di anti windup). 272 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove Si noti (in Figura 33.20) inoltre che chiudendo il parzializzatore le caratteristiche del sistema (tempo di salita, di assestamento ed overshoot) restano perlopiù invariate. Possiamo quindi concludere che il controllo progettato stabilisce un buon compromesso tra le caratteristiche analizzate (overshoot, undershoot e tempo di salita e di assestamento) risultando anche robusto a disturbi esterni, fino ad un limite dettato da cause fisiche non controllabili (si veda il caso in cui il parzializzatore è a meno di 30° ed il segnale di controllo cerca di diventare negativo). di C.Bettini, A.Guiggiani e T.Lorenzetti 273 Prove 274 di C.Bettini, A.Guiggiani e T.Lorenzetti Figura 33.5.: VI Labview del sistema di controllo Prove Figura 33.6.: Risposta Simulink con overshoot Figura 33.7.: Risposta reale con overshoot Figura 33.8.: Caratteristiche della risposta reale di C.Bettini, A.Guiggiani e T.Lorenzetti 275 Prove Figura 33.9.: Rispsosta Simulink senza overshoot Figura 33.10.: Rispsosta reale senza overshoot 276 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove Figura 33.11.: Caratteristica della reale Figura 33.12.: Risposta ottima in Simulink di C.Bettini, A.Guiggiani e T.Lorenzetti 277 Prove Figura 33.13.: Rispsosta reale Figura 33.14.: Risposta reale (particolare) 278 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove Figura 33.15.: Rispsosta con disturbo (nei cerchi sono presenti i disturbi dettati dalla chiusura del parzializzatore) Figura 33.16.: Controllo con disturbo (nei cerchi sono presenti i disturbi dettati dalla chiusura del parzializzatore) di C.Bettini, A.Guiggiani e T.Lorenzetti 279 Prove Figura 33.17.: Risposta con disturbo (nei cechi sono presenti i disturbi dettati dalla chiusura del parzializzatore) Figura 33.18.: Controllo con disturbo (nei cechi sono presenti i disturbi dettati dalla chiusura del parzializzatore) 280 di C.Bettini, A.Guiggiani e T.Lorenzetti Prove Figura 33.19.: Saturazione del segnale di controllo per valori di tensione negativi Figura 33.20.: Particolare con parzializzatore a 80° di C.Bettini, A.Guiggiani e T.Lorenzetti 281 Parte VII. Elaborati di stima e identificazione 283 34. Applicazioni del filtro di Benes 34.1. Introduzione In questo elaborato si analizzano le prestazioni di tre differenti tecniche per risolvere il problema di filtraggio di un sistema stocastico non lineare. In particolare vengono analizzati i filtri: 1. Filtro di Benes 2. Extended Kalman Filter (EKF) 3. Particle Filter Il confronto delle prestazioni viene eseguito tramite le simulazioni Monte Carlo analizzando la deviazione standard dell’errore di stima. Dato un generico sistema stocastico non lineare con rumori additivi in forma: ẋ = f (xt , ut ) + wt y = h(xt ) + vt con: x0 ⊥ (wt ⊥ vt ) wt ∼ wn (0, Qt ) 285 34.1 Introduzione vt ∼ wn (0, Rt ) È noto che il problema di filtraggio Bayesiano per tale sistema ammette soluzione sotto forma di filtro a dimensione finita nei casi: • caso lineare Gaussiano; • caso TD; • caso di Benes, in cui l’equazione di uscita h(xt ) è lineare, i disturbi wt , vt sono Gaussiani e la funzione non lineare f (xt ) soddisfa: " # af 0 0 0 + f (xt ) f (xt ) = x Ax + b x + c tr ax con ARnxn , bRn e cRn con n pari alla dimensione del vettore di stato x. Si noti che nel caso scalare la relazione diventa: df + f 2 (x) = ax2 + bx + c dx soddisfatta dalla funzione: f (x) = tanh(x) = ex − e−x ex + e−x poste a = b = 0 e c = 1. Si va così a descrivere il sistema dinamico non lineare in studio: 286 ẋ = tanh(xt ) + wt y = xt + v t di C.Bettini, A.Guiggiani e T.Lorenzetti (34.1) 34.1 Introduzione Si noti che la stima migliore di xt , nel senso del minimo errore quadratico medio (MEQM), è data da: ˆ x̂t/τ = E {xt |Yτ } = x p(xt , t | Yτ ) dx dove Yτ rappresenta le osservazioni {y0 , y1 , ..., yτ } e la funzione p(·) è la densità di probabilità di x data dalle osservazioni Yτ . di C.Bettini, A.Guiggiani e T.Lorenzetti 287 34.2 Modello del sistema 34.2. Modello del sistema Lo schema Simulink del sistema Equazione 34.1 è mostrato in Figura 34.1: Figura 34.1.: Schema Simulink del sistema dinamico non lineare Nelle figure Figura 34.2 e Figura 34.3 si analizza l’andamento dello stato e della sua derivata in assenza di rumori a partire da diverse condizioni iniziali dello stato (x0 ). Si noti inoltre come la derivata rappresenti una saturazione e il diverso comportamento del sistema alle condizioni iniziali. 288 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.3 Filtro di Benes Figura 34.2.: Andamento dello stato del sistema in assenza di rumori 34.3. Filtro di Benes La struttura del problema di Benes è data dalle equazioni: dxt dt yt = f (xt ) + wt = xt + vt con wt , vt rumori bianchi Gaussiani a media nulla e varianza rispettivamente q, r > 0. Se la funzione f (·) soddisfa la condizione: f˙ + f 2 = ax2 + bx + c considerati a, b, c costanti ed a ≥ 0, allora la stima ottima nel senso MMSE x̂t/τ = E [xt |Yτ ] esiste e dipende da un numero finito di elementi della funzione densità di probabilità di x̂t/τ . Si dimostra che almeno una di queste funzioni esiste ed è esprimibile come: di C.Bettini, A.Guiggiani e T.Lorenzetti 289 34.3 Filtro di Benes Figura 34.3.: Andamento della derivata dello stato in assenza di rumori f (x) = Aex − Be−x Aex + Be−x e che, inoltre, nel caso particolare in cui A = B, si ottiene: f (x) = tanh(x) che è la funzione in studio tipicamente utilizzata nella modellazione di saturazioni in apparaecchi elettronici. Si ottiene così un’applicazione del filtro di Benes al sistema appena descritto come in Figura 34.4. Per quanto riguarda lo studio del filtro di Benes si consideri il sistema a tempo disreto (TD). 290 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.3 Filtro di Benes Figura 34.4.: Diagramma a blocchi del filtro di Benes Per discretizzare si fa riferimento all’approssimazione di Eulero del sistema, quindi: xk+1 yk = fk (xk ) + w̄k = xk + vk (34.2) dove w̄k = Ts wk ∼ wn (0, Ts2 Q) wk , vk sono mutualmente indipendenti con distribuzioni gaussiane a media nulla e varianza rispettivamente Qk ed Rk . Per quanto riguarda la funzione non lineare, essa si traduce come: fk (xk ) = xk + Ts · f (xk ) dove f (xk ) è la funzione non lineare calcolata in xk e Ts è il tempo di campionamento adottato. di C.Bettini, A.Guiggiani e T.Lorenzetti 291 34.3 Filtro di Benes Si ha così che la stima ottima nel senso MMSE è data (in accordo a [?]) dalla relazione: x̂k/k = E {xk |Yk } = mk + Pk f (mk ) dove Yk sono le uscite fino all’istante k mentre Pk ed mk sono date dalle due equazioni ricorsive calcolate in accordo al diagramma in Figura 34.5. Figura 34.5.: Schema del filtro di Benes Quindi: 2 P¯k ¯ Pk = Pk − ¯ Pk + Rk 292 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.3 Filtro di Benes mk = m̄k + Pk (yk − m̄k ) Rk Per quanto riguarda la predizione, questa si può effettuare ponendo: dm̄(t) =0 dt dP̄ (t) =1 dt che, per tk−1 < t < tk rende: P¯k = (tk − tk−1 ) + Pk−1 = Ts + Pk−1 m̄k = mk−1 Il filtro appena descritto è inizializzato con due valori m0 , P0 che sono rispettivamente la media e la varianza dello stato iniziale avente distribuzione gaussiana: x0 v N (m0 , P0 ) Caso ideale di C.Bettini, A.Guiggiani e T.Lorenzetti 293 34.3 Filtro di Benes Si va adesso ad analizzare le prestazioni del filtro di Benes nel caso relativamente ideale, ovvero ponendo le varianze associate ai rumori q ed r, pari a: q = 0.02 r = 0.02 Si ottengono così le prestazioni riportate in Figura 34.6 e seguenti. Figura 34.6.: Stima con filtro di Benes Rumore sullo stato Si va adesso ad analizzare le prestazioni del filtro di Benes nel caso in cui si è in presenza di un accentuato disturbo di processo. Si pone in simulazione le varianze associate ai rumori q ed r, rispettivamente rumore di processo e di misura, pari a: q=1 294 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.3 Filtro di Benes Figura 34.7.: Errore quadratico medio con filtro di Benes r = 0.02 Si analizzano in Figura 34.8 e seguenti le prestazioni del filtro di Benes. Figura 34.8.: Stato e uscita nel caso di rumore sullo stato a elevata varianza di C.Bettini, A.Guiggiani e T.Lorenzetti 295 34.3 Filtro di Benes Figura 34.9.: Stima e stato reale nel caso di rumore sullo stato a elevata varianza Rumore sull’uscita Si consideri adesso la situazione opposta: un rumore di misura molto elevato ed uno sullo stato minore. Si pone così: q = 0.02 r=1 e si osservi le prestazioni del filtro in Figura 34.11 e seguenti. 296 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.4 Filtraggio EKF Figura 34.10.: Errore quadratico medio nel caso di rumore sullo stato a elevata varianza 34.4. Filtraggio EKF La logica del filtraggio EKF si basa sull’applicazione del filtro di Kalman al sistema linearizzato nel punto di interesse. Definito il vettore xk caratteristico dello stato all’istante k, si va a discretizzare il sistema. Nel caso generale considerato il sistema precedentemente esposto Equazione 34.1 si calcola: xt+1 yt = f (x̂(t|t), ut ) + = h(x̂(t|t − 1)) + df (x,ut ) dx x=x̂(t|t) dh(x) dx x=x̂(t|t−1) xet/t + o(xn ) + wt xet/t−1 + o(xn ) + vt con x̂(t|t) pari alla stima dello stato basata su y t e si definisce l’errore di predizione: xe(t|t) = x(t) − x̂(t|t) Quindi considerata la funzione: di C.Bettini, A.Guiggiani e T.Lorenzetti 297 34.4 Filtraggio EKF Figura 34.11.: Stato e uscita nel caso di rumore di misura a elevata varianza ex − e−x f (x) = tanh(x) = x e + e−x si calcola: df (x) = Ak = 1 + tanh(x)2 dx dh(x) = Ck = 1 dx Algoritmo Per quanto riguarda l’inizializzazione, si pongono le variabili (x1/0 ; P1/0 ). L’algorimo si suddivide in due fasi: correzione e predizione ( Figura 34.14). Correzione 298 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.4 Filtraggio EKF Figura 34.12.: Stima e stato reale nel caso di rumore di misura a elevata varianza Ck = ∂hk (x̂k/k−1 ) ∂x 0 Sk = Rk + Ck Pk/k−1 Ck 0 Lk = Pk/k−1 Ck Sk−1 ek = yk − hk x̂k/k−1 x̂k/k = x̂k/k−1 + Lk ek 0 Pk/k = Pk/k−1 − Lk Sk Lk Predizione Ak = ∂fk x̂ ∂x k/k−1 x̂k+1/k = fk (x̂k/k , uk ) 0 Pk+1/k = Ak Pk/k Ak + Qk dove fk (x̂k/k , uk )è la funzione discretizzata col metodo di Eulero Equazione 34.2. Caso ideale di C.Bettini, A.Guiggiani e T.Lorenzetti 299 34.4 Filtraggio EKF Figura 34.13.: Errore quadratico medio nel caso di rumore di misura a elevata varianza Figura 34.14.: Passi caratterizzanti il filtraggio EKF Si va adesso ad analizzare le prestazioni dell’EKF nel caso relativamente ideale, ovvero ponendo le varianze associate ai rumori q ed r, pari a: q = 0.02 r = 0.02 Si ottengono così le prestazioni riportate in Figura 34.15 e seguenti. 300 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.4 Filtraggio EKF Figura 34.15.: Stima nel caso di filtraggio EKF Rumore sullo stato Come nel caso del filtraggio di Benes, analizziamo il problema con rumori di processo molto incisivi. Posti q = 1 e r = 0.02 osserviamo in Figura 34.17 e successive le prestazioni del filtraggio EKF. Rumore sull’uscita Si osservino infine in Figura 34.20 e successive le prestazioni del filtro EKF nel caso opposto: rumore di processo poco intenso (q = 0.02) e rumore di misura molto incisivo (r = 1). di C.Bettini, A.Guiggiani e T.Lorenzetti 301 34.5 Filtraggio PF Figura 34.16.: Errore quadratico medio nel caso di filtraggio EKF 34.5. Filtraggio PF Il filtro PF (particle filter) equivale ad una implementazione campionata del filtro Bayesiano sfruttando il metodo Monte Carlo per l’integrazione (IMC ). Invece di una soluzione analitica o di una approssimazione numerica, infatti, il filtro PF considera una grossa mole di dati al fine di stimare al meglio la densità di probabilità (ddp) dello stato. L’idea che sta alla base del metodo è quella di rappresentare la ddp con una serie di campioni casuali (da qui il nome particles) tali che, incrementandone il numero, sia possibile ottenere una rappresentazione stimata della ddp desiderata. Algoritmo Inizializzazione: Si inizializzano i campioni di xi e i pesi w̄i come : xi (0) ∼ N (x(0); P (0)) 302 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.5 Filtraggio PF Figura 34.17.: Stato e uscita nel caso di rumore sullo stato a elevata varianza w̄i (0) = 1 M calcolati rispettivamente per i = 1, ... , M , dove M è il numero di particelle scelte dal progettista. Esecuzione: per t > 1..... e per i = 1, ... , M si calcolano i campioni: xi (t) v q(x(t)|xi (t − 1), y(t)) ed i pesi: wi (t) = wi (t − 1) p(y(t)|xi (t)) p(xi (t)|xi (t − 1)) q(xi (t)|xi (t − 1), y(t)) Si vanno a normalizzare i pesi considerando: di C.Bettini, A.Guiggiani e T.Lorenzetti 303 34.5 Filtraggio PF Figura 34.18.: Stima e stato nel caso di rumore sullo stato a elevata varianza wi (t) w̄i (t) = PM j=1 wj (t) Si va a calcolare la stima dello stato: x̂(t) = M X w̄i (t) xi (t) i=1 Si stima la covarianza come: P̂ (t) = M X 0 w̄i (t) [xi (t) − x̂(t)] [xi (t) − x̂(t)] i=1 Si noti che nell’algoritmo il progettista deve scegliere la funzione di importanza q(x(t)|x(t − 1), y(t)), la scelta più semplice è porre: q(x(t)|x(t − 1), y(t)) = p(xi (t)|xi (t − 1)) 304 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.5 Filtraggio PF Figura 34.19.: Errore quadratico medio nel caso di rumore sullo stato a elevata varianza L’unico vincolo sulla funzione di importanza è che il suo supporto contenga quello di p(x(t)|y t ); si sceglie così q(·) come una distribuzione Gassiana. Ricampionamento Il fenomeno della degenerazione dei pesi, che è causato dalla fattorizzazione della funzione di importanza, è quantificabile valutando il numero dei campioni effettivi Nef f i cui pesi sono significativi: 1 Nef f = PN i=1 w̄i2 La degenerazione degrada l’accuratezza della stima e non utilizza in modo efficiente le risorse di calcolo. Si devono, infatti, sprecare risorse di calcolo per propagare campioni con pesi nulli che, non contribuiscono alla stima. Tale fenomeno può essere contrastato effettuando un ricampionamento. Con un ricampionamento si vanno a re-impostare i valori dei pesi con una certa logica. Nel nostro caso si è scelto di porre, se vale la condizione di ricampionamento, i pesi: di C.Bettini, A.Guiggiani e T.Lorenzetti 305 34.5 Filtraggio PF Figura 34.20.: Stato e uscita nel caso di rumore di misura a elevata varianza wi (t) ∼ N (0, q) ovvero distribuendo i pesi in modo Gaussiano con varianza q. Caso ideale Si va adesso ad analizzare le prestazioni del PF nel caso relativamente ideale, ovvero ponendo le varianze associate ai rumori q ed r, pari a: q = 0.02 r = 0.02 Si ottengono così i risultati mostrati in Figura 34.23 e seguenti. Si noti in Figura 34.25 l’andamento dei campioni dell’algoritmo PF in uso. 306 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.5 Filtraggio PF Figura 34.21.: Stima e stato reale nel caso di rumore di misura a elevata varianza Rumore sullo stato Anche in quest’ultimo caso si analizzano prima le prestazioni del filtro con rumore di processo molto elevato (q = 1, r = 0.02), poi nel caso di rumore di misura incisivo (q = 0.02, r = 1). Essendo il filtraggio tramite filtro particellare caratterizzato anche dal numero di campioni utilizzati, si riportano i risultati ( Figura 34.26 e successive) per N = 10, N = 250 ed N = 5000. Si noti infine che per N = 10 è necessario applicare il ricampionamento (evidenziato in rosso) a causa della degradazione della stima, osservando come si comporta l’errore, sopratutto dopo il primo ricampionamento. di C.Bettini, A.Guiggiani e T.Lorenzetti 307 34.5 Filtraggio PF Figura 34.22.: Errore quadratico medio nel caso di rumore di misura a elevata varianza Figura 34.23.: Stima nel caso di filtraggio PF 308 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.5 Filtraggio PF Figura 34.24.: Errore quadratico medio nel caso di filtraggio PF Figura 34.25.: Evoluzione dei campioni nell’algoritmo PF, caso N = 250 di C.Bettini, A.Guiggiani e T.Lorenzetti 309 34.5 Filtraggio PF Figura 34.26.: Stato e uscita nel caso di rumore sullo stato a elevata varianza N = 10 Figura 34.27.: Stima e stato reale nel caso di rumore sullo stato a elevata varianza N = 10 310 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.5 Filtraggio PF Figura 34.28.: Errore quadratico medio nel caso di rumore sullo stato a elevata varianza N = 10 Figura 34.29.: Stima e stato reale nel caso di rumore sullo stato a elevata varianza N = 250 di C.Bettini, A.Guiggiani e T.Lorenzetti 311 34.5 Filtraggio PF Figura 34.30.: Errore quadratico medio nel caso di rumore sullo stato a elevata varianza N = 250 Figura 34.31.: Stima e stato reale nel caso di rumore sullo stato a elevata varianza N = 5000 312 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.5 Filtraggio PF Figura 34.32.: Errore quadratico medio nel caso di rumore sullo stato a elevata varianza N = 5000 di C.Bettini, A.Guiggiani e T.Lorenzetti 313 34.5 Filtraggio PF Rumore sull’uscita Si riporta in Figura 34.33 e successive l’ultimo caso, in cui il rumore sull’uscita (caratterizzato da una varianza r = 1) è più incisivo rispetto a quello di processo (q = 0.02). Figura 34.33.: Stato e uscita nel caso di rumore di misura a elevata varianza N = 10 314 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.5 Filtraggio PF Figura 34.34.: Stima e stato reale nel caso di rumore di misura a elevata varianza N = 10 Figura 34.35.: Errore quadratico medio nel caso di rumore di misura a elevata varianza N = 10 di C.Bettini, A.Guiggiani e T.Lorenzetti 315 34.5 Filtraggio PF Figura 34.36.: Stima e stato reale nel caso di rumore di misura a elevata varianza N = 250 Figura 34.37.: Errore quadratico medio nel caso di rumore di misura a elevata varianza N = 250 316 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.5 Filtraggio PF Figura 34.38.: Stima e stato reale nel caso di rumore di misura a elevata varianza N = 5000 Figura 34.39.: Errore quadratico medio nel caso di rumore di misura a elevata varianza N = 5000 di C.Bettini, A.Guiggiani e T.Lorenzetti 317 34.6 Prove Monte Carlo 34.6. Prove Monte Carlo Le prestazioni dei filtri presentati sono valutate attraverso l’analisi dell’errore quadratico medio di stima con il metodo Monte Carlo. Il metodo Monte Carlo si basa sull’effettuare M prove indipendenti tra loro in termini di condizione iniziale e rumore. L’indipendenza del rumore tra una prova e l’altra è garantita dall’aver impostato randseed come nel campo Initial seed del blocco Simulink che genera i rumori gaussiani. Date M prove la deviazione standard (σet ) dell’errore di stima è così definita: σet = v u u t M 1 X (xt − x̂t/t )2 M t=1 Posto M = 500, si riporta di seguito le prestazioni dei filtri in termini di deviazione standard dell’errore di stima calcolate attraverso il metodo Monte Carlo al variare della varianza q (caratteristica del rumore di processo) e della varianza r (relativa al rumore di misura). Nelle simulazioni il filtro PF è considerato con un numero di campioni N pari a 250. Le simulazioni di durata 30s sono inizializzate con uno stato iniziale dei filtri definito come: x0 ∼ N (15; 5) mentre lo stato iniziale vero è impostato pari a x0 = 5 Caso ideale q = 0.2 318 r = 0.2 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.6 Prove Monte Carlo Figura 34.40.: EQM - caso ideale Figura 34.41.: EQM - caso ideale Rumore sullo stato q = 1 r = 0.2 Figura 34.42.: EQM - rumore sullo stato di C.Bettini, A.Guiggiani e T.Lorenzetti 319 34.6 Prove Monte Carlo Figura 34.43.: EQM - rumore sullo stato Rumore sull’uscita q = 0.2 r=1 Figura 34.44.: EQM - rumore sull’uscita Figura 34.45.: EQM - rumore sull’uscita 320 di C.Bettini, A.Guiggiani e T.Lorenzetti 34.6 Prove Monte Carlo Come si può notare dalla lettura delle figure, il filtro di Benes ha ottime prestazioni in tutte le casistiche prese in considerazione, da notare che nel caso di forti disturbi sull’uscita si ha un EQM più alto. Il PF normalmente ha un comportamento peggiore rispetto agli altri filtri, ma nel caso di forti rumori sull’uscita, a fronte di una maggiore onerosità computazionale, si comporta meglio dell’EKF (il quale si comporta peggio nel caso di forti disturbi in uscita). Inoltre è evidente come il PF ha una convergenza più lenta rispetto agli altri filtri. Si riporta di seguito altre simulazioni considerando un numero di campioni del PF pari a 2000. Figura 34.46.: EQM - caso ideale Figura 34.47.: EQM - rumore sullo stato di C.Bettini, A.Guiggiani e T.Lorenzetti 321 34.6 Prove Monte Carlo Figura 34.48.: EQM - rumore sull’uscita 322 di C.Bettini, A.Guiggiani e T.Lorenzetti 35. Stima dello stato in presenza di osservazioni binarie L’elaborato si incentra sull’implementazione di algoritmi di filtraggio ricorsivo applicati ad un sistema fisico reale non lineare, per stimare lo stato dell’impianto non direttamente accessibile. Il processo oggetto dello studio è composto dalla cascata di due vasche di differenti dimensioni e portate. La portata di ingresso della prima vasca è regolata da un segnale non costante generato in modo che il livello d’acqua di entrambe non superi mai un valore massimo prefissato. L’ingresso della seconda invece è determinato dalla portata di uscita della prima. L’impianto reale è stato modellato tramite l’utilizzo di Simulink e discretizzato in modo tale da poter implementare in Matlab gli algoritmi EKF ed UKF necessari alla stima del suo stato. Il primo si basa sul Filtro di Kalman Esteso mentre il secondo si basa sulla Trasformazione Unscented. L’impianto presenta in uscita un sensore che fornisce i dati necessari all’elaborazione della stima. Inizialmente ‘e stato implementato un sensore a due livelli; il problema ‘e stato poi generalizzato utilizzando un sensore con un numero finito di livelli. In conclusione dell’elaborato ‘e stato effettuato il confronto tra gli algoritmi utilizzati tramite la valutazione dell’errore quadratico medio basata su analisi MonteCarlo. 35.1. Modellazione Il sistema è descritto dalle equazioni non lineari: 323 35.1 Modellazione √ = Ku − A1 2gx1 + w1 √ √ S2 ẋ2 = A1 2gx1 − A2 2gx2 + w2 S1 ẋ1 le quali descrivono la relazione tra l’ingresso u di alimentazione ed i livelli dell’acqua nelle due vasche. E’ stato necessario ricavare uno schema Simulink ( Figura 35.1) in grado di simulare il sistema dinamico reale tenendo conto di un errore di processo riguardante lo stato ed un errore di misura riguardante l’uscita. L’ingresso scelto è di tipo sinusoidale ed è stato ricavato a seguito di una serie di test sul sistema dinamico simulato in modo tale che gli stati x1 ed x2 siano compresi tra un livello minimo nullo ed un livello massimo unitario h e che l’altezza della seconda vasca oscilli intorno al valor medio h2 . Figura 35.1.: Modello Simulink delle due vasche I parametri del sistema dinamico sono: 324 di C.Bettini, A.Guiggiani e T.Lorenzetti 35.1 Modellazione S1 = 0.5 = Area di base della prima vasca; S2 = 0.25 = Area di base della seconda vasca; A1 = 0.019 = Superficie della prima sezione di scarico; A2 = 0.016 = Superficie della seconda sezione di scarico; K = 0.035 = Costante di proporzionalità caratteristica della pompa; g = 9.81 = Accelerazione di gravità; h = 1 = Altezza massima del livello dell’acqua; v = 0.01 = deviazione standard del rumore di misura; w1 , w2 = 0 = errore di processo. L’uscita y del sistema è definita come lo stato della seconda vasca a cui si somma un rumore di misura v: y = x2 + v Per quanto riguarda il modello del sensore lineare a due livelli, esso fornisce tramite un comparatore misure binarie bk pari a: bk = 0, se y < h 2 1, se y > h 2 Successivamente è stato associato al valore booleano bk un valore yk così definito: yk = h h, 0 1 xk + vk = 4 3h i 4 se bk = 0 , se bk = 1 dove vk rappresenta il rumore additivo introdotto dal comparatore: δ2 var (vk ) = σ + 3 2 di C.Bettini, A.Guiggiani e T.Lorenzetti 325 35.2 Algoritmo EKF dove δ = 2lh , con h pari all’altezza massima ed l pari al numero di livelli del comparatore. 35.2. Algoritmo EKF Il primo algoritmo di filtraggio non lineare utilizzato per la stima dello stato del sistema dinamico è EKF. Formalmente l’algoritmo EKF è così definito: per k = 1, 2, ... Correzione: Ck = ∂hx ∂x x=x̂ k|k−1 0 Sk = Rk + Ck Pk|k−1 Ck 0 Lk = Pk|k−1 Ck Sk−1 x̂k|k = x̂k|k−1 + Lk ek 0 Pk|k = Pk|k−1 − Lk Sk Lk Predizione: Ak = ∂fk ∂x x=x̂ k|k−1 Dk = gk x̂k|k x̂k+1|k = fk x̂k|k , uk 0 0 Pk+1|k = Ak Pk|k Ak + Dk Qk Dk Per implementare l’algoritmo è necessaria la discretizzazione del sistema e la ricerca di una matrice Ak definita come il giacobiano della dinamica discretizzata dello stato assumendo un tempo di campionamento δk = 0.1. Definito il vettore di stato: xk = 326 h xk (1) xk (2) i0 di C.Bettini, A.Guiggiani e T.Lorenzetti 35.2 Algoritmo EKF tramite la discretizzazione di Eulero si ottiene: √ 2gxk (1) x (1) + δ Ku − A k k 1 S1 √ √ xk+1 = 2gxk (1) 2gxk (2) − A2 S2 xk (2) + δk A1 S1 Discretizzando è opportuno considerare un disturbo di processo (dovuto all’approssimazione introdotta dal metodo) di cui si tiene conto attraverso: 10−4 0 Q= 0 10−4 Si può quindi ricavare la matrice Ak : √ 1 Ak = 2gxk (1) + δk √ S1 A g 2gx (1) δk 1 S1 k A1 g 0 1 − δk A2 g √ 2gxk (2) S2 Le condizioni iniziali per questo specifico sistema sono definite in modo tale da avere uno stato di partenza prossimo a zero ed una varianza molto elevata data la scarsa quantità di informazioni inerenti lo stato. A tale proposito si ‘e posto: x1|0 = P1|0 h 0.01 0.01 i0 0.5 0 = 0 0.5 Data la linearità del sensore in uso, la matrice Ck definita nell’algoritmo non va ricavata poichè costante e più precisamente pari a: Ck = h 0 1 i L’algoritmo EKF è stato implementato in tre programmi diversi che si differenziano nel modo in cui viene trattato il segnale di ingresso all’algoritmo. di C.Bettini, A.Guiggiani e T.Lorenzetti 327 35.2 Algoritmo EKF Nel primo programma viene utilizzato un comparatore a due livelli implementato direttamente in Simulink. L’associazione tra misura binaria bk uscente dal comparatore e misura continua yk è effettuata a livello Matlab. Il secondo programma implementa un comparatore a l-livelli al cui ingresso è posto il segnale Simulink y a cui si associa in uscita il segnale yk . Il terzo, implementato anch’esso completamente in Matlab, realizza il caso l = ∞ ponendo in ingresso all’algoritmo direttamente l’uscita y del sistema non comparata. Dalle simulazioni effettuate si ottengono i risultati mostrati in Figura 35.2 e seguenti. Figura 35.2.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a due livelli Figura 35.3.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a quattro livelli 328 di C.Bettini, A.Guiggiani e T.Lorenzetti 35.2 Algoritmo EKF Figura 35.4.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a otto livelli Figura 35.5.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a sedici livelli di C.Bettini, A.Guiggiani e T.Lorenzetti 329 35.2 Algoritmo EKF Figura 35.6.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a 2048 livelli Figura 35.7.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a infiniti livelli 330 di C.Bettini, A.Guiggiani e T.Lorenzetti 35.3 Algoritmo UKF 35.3. Algoritmo UKF Il secondo algoritmo di filtraggio non lineare utilizzato è UKF, basato sulla Trasformazione Unscented, che inizializza dei campioni aggiornati in modo ricorsivo al fine di stimare lo stato del sistema. Inizialmente viene definito il parametro na = nx + nw + nv pari alla somma delle dimensioni di stato, sistema ed errori (di processo e di misura). Nel caso in esame, essendo i rumori additivi e gaussiani, na è pari alla dimensione dello stato (quindi na = 2). Definiamo le condizioni iniziali ponendo lo stato nullo e la varianza elevata: x1|0 = h 0 0 i0 P1|0 0.5 0 = 0 0.5 Inizializziamo poi i parametri necessari all’algoritmo: w0 = 1 − n30 , pari al peso del q nx 0 campione centrale; wi = 1−w pari al peso degli altri campioni; s = . 2nx 1−w0 Formalmente l’algoritmo UKF è così definito: per k = 1, 2, ... x̂a = P̂ a = x̂t−1|t−1 0nw 0nv Pt−1|t−1 0 0 0 Q 0 0 0 R di C.Bettini, A.Guiggiani e T.Lorenzetti 331 35.3 Algoritmo UKF dove Q, R sono le varianze dei rumori rispettivamente di processo e di misura. Nel nostro caso particolare, essendo i rumori additivi e gaussiani, si ha che: x̂a = x̂t−1|t−1 P̂ a = Pt−1|t−1 Calcoliamo così i campioni: χa0 = x̂a χa1 = x̂a + sTia χa2 = x̂a − sTia χa = h χ x χ w χ v i0 = χx nel nostro caso Predizione: χxi = ft (χxi , χw i , ut ) x̂t|t−1 = n X wi χxi i=0 Pt|t−1 = n X wi x̂t|t−1 − χxi x̂t|t−1 − χxi 0 i=0 ζi = ht (χxi , χvi ) ŷt|t−1 = n X wi ζi i=0 Correzione: 332 di C.Bettini, A.Guiggiani e T.Lorenzetti 35.3 Algoritmo UKF St = Pn i=0 wi ŷt|t−1 − ζi Ptxy = Pn x i=0 wi x̂t|t−1 − χi ŷt|t−1 − ζi 0 ŷt|t−1 − ζi 0 Lt = Ptxy St−1 x̂t|t = x̂t|t−1 + Lt yt − ŷt|t−1 0 Pt|t = Pt|t−1 − Lt St Lt Analogamente al caso EKF l’algoritmo è stato implementato in Matlab sfruttando il modello discretizzato del sistema dinamico. Le simulazioni sono state effettuate in primo luogo sfruttando un sensore a due livelli ed in secondo luogo generalizzando il problema fino ad arrivare al caso di infiniti livelli. I risultati ottenuti sono riportati in Figura 35.8 e seguenti. Figura 35.8.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a due livelli di C.Bettini, A.Guiggiani e T.Lorenzetti 333 35.3 Algoritmo UKF Figura 35.9.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a quattro livelli Figura 35.10.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a otto livelli 334 di C.Bettini, A.Guiggiani e T.Lorenzetti 35.3 Algoritmo UKF Figura 35.11.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a sedici livelli Figura 35.12.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a 2048 livelli di C.Bettini, A.Guiggiani e T.Lorenzetti 335 35.3 Algoritmo UKF Figura 35.13.: Stati x1 ed x2 reali e loro stima filtrata nel caso di comparatore a infiniti livelli 336 di C.Bettini, A.Guiggiani e T.Lorenzetti 35.4 Conclusioni 35.4. Conclusioni Le prestazioni dei due algoritmi sono state valutate attraverso l’analisi dell’MSE con il metodo MonteCarlo. Date M prove, si indica con i la i-esima realizzazione e si ricava l’MSE per ogni istante temporale t come: M SEt = M 0 1 X xi − x̂it|t xi − x̂it|t M i=0 Per implementare questo metodo di analisi è stato scelto di eseguire un’unica simulazione del sistema dinamico prelevando lo stato vero x2 a cui è stato aggiunto, a livello Matlab, un rumore gaussiano a media nulla (con v = 10−2 ) generato in modo indipendente per ognuna delle 50 realizzazioni. Si è così realizzato un quarto file Matlab grazie al quale è possibile analizzare l’evoluzione temporale dell’MSE dei due algoritmi (rispettivamente per 2, 4, 8 e 16 livelli) come mostrato in Figura 35.14 e Figura 35.15. Dai risultati ottenuti si evince che l’MSE diminuisce sensibilmente all’aumentare dei livelli per entrambi gli algoritmi ma in modo più marcato nel caso EKF. Le prestazioni dei due algoritmi sono state inoltre confrontate attraverso l’analisi temporale dell’errore quadratico medio definito come: M SE = N 0 1 X xt − x̂t|t xt − x̂t|t N t=0 I risultati ottenuti, per ogni algoritmo da l = 2 ad l = ∞, sono riportati nella tabella Tabella 35.1. Anche in questo caso si può notare come l’algoritmo EKF sia migliore rispetto a quello UKF e come per entrambi l’MSE diminuisca all’aumentare dei livelli. di C.Bettini, A.Guiggiani e T.Lorenzetti 337 35.4 Conclusioni Figura 35.14.: Andamento dell’errore quadratico medio calcolato tramite analisi MonteCarlo; algoritmo EKF, comparatore a 2 (a), 4 (b), 8 (c) e 16 (d) livelli Tabella 35.1.: Errore Quadratico Medio temporale calcolato nei due metodi utilizzati 338 di C.Bettini, A.Guiggiani e T.Lorenzetti 35.4 Conclusioni Figura 35.15.: Andamento dell’errore quadratico medio calcolato tramite analisi MonteCarlo; algoritmo UKF, comparatore a 2 (a), 4 (b), 8 (c) e 16 (d) livelli di C.Bettini, A.Guiggiani e T.Lorenzetti 339