MATLAB-SIMULINK€¦ · 3 Finestra di avvio (v. 7.8.0) Cartella corrente Avvio SIMULINK Editor...

MATLAB-SIMULINK

Simulink

Ing. Alessandro Pisano

pisano@diee.unica.it

IntroduzioneLibrerie e blocchi elementari (1) Realizzazione di un modello Esempio: costruzione e visualizzazione di una sinusoide Scelta del solutore Modifica delle impostazioni predefinite Utilizzo di variabili dal workspaceLibrerie e blocchi elementari (2) Esportazione dati verso il Workspace e su file esterno Esempio: filtro passa basso Integratore Esecuzione automatizzata di test Analisi spettrali (FFT) Filtraggio digitale Esempio: sistema termico ad 1 e 2 gradi di libertàUtilizzo di blocchi Trasnfer functionEsempio: Simulazione di un sistema di regolazione di temperaturaCreazione di sottosistemiMaskEsempio Distribuzione di temperatura nel rotore di una turbinaRealizzazione di sistemi MIMO LTIVariabili popup e checkboxEsempio: sistema di frenatura con ABSEmbedded Matlab FunctionToolbox avanzati. SimMechanics e SimDriveline

379111620272831374043445963687381859099109112123133

Indice

Alessandro Pisano - pisano@diee.unica.it

Finestra di avvio (v. 7.8.0)

Cartella corrente

Avvio SIMULINK

Editor M-files

Programmazione dei modelli di simulazione per via grafica

Possibilità di suddividere il modello complessivo in sottosistemi paralleli o embricati (un sottosistema può essere importato direttamente in un modello di simulazione differente)

Possibilità di definire finestre di parametrizzazione (Masks)

Esportazione nel workspace Matlab dei risultati della simulazione

Esecuzione automatizzata di test

Toolbox avanzati: SimDriveline e SimMechanics

Avvio SIMULINK

Librerie principali

Libreria “Commonly used blocks”

Libreria “Sinks”

Libreria “Sources”

Pagina di lavoro

Realizzazione di un modello Simulink

1. Importare nella pagina di lavoro i blocchi elementari Simulinknecessari, trascinandoli con il mouse dalla rispettiva libreria (drag-and-drop)

3. Collegare tra loro i blocchi Simulink tracciando le opportune linee di interconnessione in modo da realizzare le funzionalità desiderate

2. Parametrizzare i blocchi Simulink nelle rispettive finestre di parametrizzazione, alle quali si accede dalla pagina di lavoro facendo doppio click con il mouse sopra il blocco stesso.

3 fasi

Esempio introduttivo: costruzione e visualizzazione di un segnale sinusoidale

Sono sufficienti due blocchi elementari: un blocco che generi il segnale desiderato, ed un blocco che ne permetta la visualizzazione.

Il primo blocco lo troveremo nella libreria “Sources” (blocco Sine Wave),Il secondo blocco (blocco Scope), si trova nella libreria “Sinks”.

I blocchi necessari vannoimportati nella pagina di lavoroUntitled trascinando con ilmouse (drag-and-drop) l’iconadel blocco all’interno dellapagina di lavoro. Il risultato ditale procedura è mostrato inFigura.

Salvare il modello e attribuire unnome al file con estensione .mdl

Si deve ora collegare l’uscita del generatore di funzione “Sine Wave” conl’ingresso del blocco di visualizzazione “Scope”.

Per effettuare un collegamento tra due blocchi vi è una procedura rapida. Sideve selezionare il blocco di origine (cliccandovi sopra), e si devesuccessivamente selezionare il blocco di destinazione con il tasto ctrlpremuto.

Un collegamento correttamente eseguito vieneindicato come in Figura

In alternativa, si può portare la freccia del mouse nel punto di origine del collegamento e quindi “tracciarlo” tenendo premuto il tasto sinistro del mouse, portandosi fino al punto di destinazione.

Devono ora essere impostati i parametri di ampiezza, frequenza e sfasamentoche definiscono la particolare sinusoide che si desidera generare. A tal fine ènecessario fare doppio click sul blocco “Sine Wave”, e come risultato si apreuna finestra di dialogo all’interno della quale vanno impostati i parametri difunzionamento.

Ampiezza

Frequenza

Sfasamento

Tasto OKAlessandro Pisano - pisano@diee.unica.it

Si deve ora impostare la durata (cioè l’intervallo temporale) della simulazione.

Cliccando sul tasto RUN vieneeseguita la simulazione.Dopo che è stata eseguita lasimulazione si può visualizzare ilsegnale generato cliccando sul bloccoScope.

La durata si può impostare direttamente dai menù della pagina di lavoro

Durata(valore di default 10.0)

Tasto RUN

Ora si aumenti la frequenza della sinusoide da 1 rad/s a 2 rad/s

Si ripeta la simulazione. Si ri-aggiorni il grafico della finestra grafica Scope

cliccando sul pulsante nella barra dei menu della finestra Scope

Grafico “spigoloso”

Si deve andare a modificare il “metododi integrazione”, che definisce il passo didiscretizzazione temporale che vieneimpiegato nella simulazione delmodello.

Il metodo di integrazione (Solver) si imposta selezionando il menu della pagina diconfigurazione Simulation->Configuration Parameters

La scelta del solutore con il quale di risolvono numericamente le equazionidifferenziali del modello è ovviamente irrilevante per il semplice esempio in esameche non coinvolge alcun legame differenziale.

Il motivo della spigolosità del grafico sta nel fatto che il metodo proposto di default(ode45 a passo variabile) ha “scelto” dei passi di discretizzazione temporalepiuttosto elevati, e sono stati quindi generati “pochi campioni” del segnale

In modelli di simulazione complessi, la scelta del solutore numerico (Runge-Kutta,Dormand-Prince, Eulero,…) e delle relative caratteristiche (passo fisso/variabile,etc) va fatta con criterio.

Runge-Kutta e Dormand-Prince sono degli ottimi solutori “general-purpose”.

Per simulazioni “multi-domain” con la copresenza di costanti di tempo moltodifferenti tra loro sono consigliati i metodi a passo variabile dedicati ai problemiStiff (es. ode15s/stiff).

Per simulazioni con elementi discontinui (non-smooth dynamics) i metodi a passovariabile talvolta forniscono risposte non veritiere. Il solutore Eulero a passo fisso,con un passo sufficientemente piccolo, è ritenuto affidabile per sistemi non-smooth.

Si scelga il solutore ode1 (Eulero) a passo fisso, e se ne imposti il Fixed-Step size a0.001

Una scelta ottimale per il solutore bilancia, per il problema in esame, la precisione

della soluzione e la mole di calcoli richiesta, che influenza il tempo di simulazione.

Per identificare un solutore adeguato serve esperienza. Il passo di campionamento

deve essere commisurato alla rapidità di variazione dei segnali in gioco. Quando si

sceglie un solutore a passo variabile si può pensare di introdurre un limite massimo

per il passo adattativo.

Quando il modello non contiene stati continui (non vi sono cioe blocchi “dinamici”

come Integratori, blocchi Transfer Fcn, etc.) Simulink usa il solutore “discrete”

anche se viene specificato un solutore differente,.

Considerazioni aggiuntive

Modifica delle impostazione predefinite per i files Simulink all’apertura

Modificare le Solver Options come in figura

Modificare le Data Import/Export Options come in figura

Si ripeta la simulazione e si riaggiorni il grafico

La sinusoide viene però mostrata a partiredall’istante t=5.Sono stati “persi” i campioni precedenti.

Il motivo è che, al fine di non saturare rapidamente la memoria del programma,vige in Simulink una impostazione di default in base alla quale nei blocchi di tipoScope vengono visualizzati e mantenuti in memoria solo gli ultimi 5000 campionidel segnale.

Il grafico della sinusoide è ora correttamenterappresentato.In base alla scelta fatta per il passo fisso delsolutore, vengono ora generati, e interpolati dalgrafico, 1000 campioni per ogni secondo dievoluzione del segnale.

Per modificare tale impostazione per uno specificoblocco si deve cliccare sul pulsante Parameters nellafinestra del blocco Scope

La finestra “Scope Parameters” ha due sottomenu: “General” e “Data History”.

Dal sottomenu Data History si deve disselezionare la check-box “Limit data points to last ..”

Ora ripetendo la simulazione e riaggiornando il grafico la sinusoide vienevisualizzata per intero

Per visualizzare un segnale costituito dalla somma di tre sinusoidi importiamonella pagina di lavoro due nuove istanze del blocco elementare Sine Wave, edimportiamo anche un blocco che rappresenti un nodo sommatore (blocco Sum

dalla libreria dei Commonly Used Blocks)

Dopo aver cancellato la linea di collegamentopreesistente tra il primo blocco Sine Wave edil blocco Scope, si realizzi la connessioneriportata in Figura.

Il blocco Sum deve esserepreliminarmente parametrizzatospecificando il numero di segnali iningresso, ed il segno con il qualeconcorrono alla sommatoria, permezzo di una stringa (es. ++−+ −+)

Scegliamo +++L’aspetto del blocco diventa

Ora si possono assegnare i parametri di ampiezza, frequenza, bias e sfasamento, delletre sinusoidi, rieseguire la simulazione e visualizzare il grafico prodotto dal bloccoScope .

I parametri di un modello possono essere espressi utilizzando delle variabili (es.A1,f1,b1, …) alle quali si può assegnare un valore con un file script Matlab daeseguirsi prima della esecuzione della simulazione.

Tutte le variabili definite nel workspace di Matlab sono disponibili e accessibili daparte dei blocchi Simulink.

Il modello può essere in questomodo riparametrizzato con estremafacilità.

phi1=0;

f2=pi;

phi2=pi/2;

f3=4*pi;

phi3=0;

Libreria “Continuous”

Libreria “Discontinuities”

Libreria “Discrete”

Libreria “Lookup Tables”

Libreria “Math operations”

Libreria “Model Verification”

Libreria “Signal routing”

Il blocco To Workspace riceve in ingresso il segnale(scalare o vettoriale) che salva nel workspace. Ilblocco si interconnette agli altri come in Figura.Per tracciare un collegamento a partire da uncollegamento preesistente si deve portare il mousenel punto di diramazione, premere il tasto destro, epoi allontanarsi e tracciare il collegamento tenendoil tasto destro premuto, fino a giungere al punto didestinazione del collegamento

Vediamo come esportare in Matlab i dati prodotti eseguendo i modelliSimulink.

Serve il blocco To Workspace dalla libreria Sinks

Il blocco To Workspace va parametrizzato specificando ilnome della variabile che verrà creata nel workspace diMatlab (scegliere y) ed il formato di salvataggio (èopportuno modificare il formato di default Structure eselezionare invece Array) .Per ottimizzare l’impiego della memoria del programma sipuò anche impostare un fattore intero di decimazione(es. con decimation = 10 i dati vengono salvati nelworkspace con uno step temporale 10 volte superiore,quindi si avranno meno elementi nel vettore y).Si mantenga il valore unitario di default.

Il vettore dei tempi viene salvato didefault sotto forma di array con ilnome tout.Bisogna però disabilitare unaimpostazione che limita a 1000 ilnumero massimo di elementi pertout. Si deve andare ne Simulation-

>Configuration Parameters, e nel

menu Data Import/Exportdisselezionare la check-box Limit

Data Points to Last

Lanciando la simulazione vengono creati nel workspace di Matlab gli array y e

Verificarlo digitando il comando whos

Si può visualizzare in Matlab il grafico del segnale con il comando

plot(tout,y),grid

Se si desidera che i dati siano non solo esportatinel workspace di Matlab ma anche salvati inmodo permanente su un file, si può generare unfile dati nel formato binario .mat con il blocco To

File, dalla libreria Sinks

I files con estensione mat sono detti mat-files.

Esportazione su file dei dati prodotti da modelli Simulink.

Il blocco deve essere parametrizzato specificandonella apposita finestra di configurazione:

-Il nome (comprensivo dell’estensione) del mat-fileche verrà creato (es. y_test1.mat).

-Il nome che verrà assegnato alla variabilequando il mat-file sarà successivamente

aperto in Matlab (scegliere y)E’ possibile impostare una decimazione dei dati.Se si sceglie una “Decimation” > 1 è bene generareanche un vettore dei tempi “sincrono” con lavariabile sottocampionata

Lanciare la simulazione, e verificare come nella cartella di lavoro sia ora presente ilfile y_test1.mat

Possono essere esportati segnali vettoriali. Per mezzo del blocco Mux (libreriaCommonly Used Blocks) si possono “aggregare” i tre segnali sinusoidali in un unicosegnale vettoriale con tre componenti.

La struttura interna dei mat-files prevede la memorizzazione dei dati in una strutturarettangolare

( ) ( ) ( ) ( )

2322212

1312111

tytytyty

I mat-files possono essere aperti successivamente in Matlab con il comando load.

>> load y_test1

Viene generata nel workspace di Matlab una variabile matriciale avente il nomespecificato nella finestra di configurazione del blocco To File e la struttura rettangolareriportata sopra.

Integratore

Nella finestra di parametrizzazione dell’integratore il parametro piu importante da settare èla condizione iniziale (Initial Condition), che di default viene impostata pari a zero

Possono essere introdotte saturazioni inferiorie/o superiori sulle uscite dell’integratore.

Se l’integratore riceve in ingresso un segnalevettoriale, genera in uscita un vettore di paridimensione che contiene l’integrale delle diversecomponenti del vettore di ingresso

Può essere anche applicato un reset sull’uscitadel’integratore. Il reset riporta l’uscitadell’integratore al valore della condizione iniziale

Integratore con ingresso vettoriale

Condizioni iniziali diverse per le uscite

Integratore con saturazione superiore

( ) ( ) ( ) ( ) ( )* sin MgtFtxktxbtxM −=++

Sistema massa molla smorzatore su piano inclinato

( ) ( )ttF 2sin1010+=

Sistema massa molla smorzatore su piano inclinato

( ) ( )ttF 2sin1010+=

( ) ( ) ( ) ( ) ( )* sin1

btx −+−−=

Eq. in forma esplicita

Modello Simulink

Modello Simulink compatto

( ) ( ) gtmtFzzktztm )()( −=+

( )temmtm 1.0

10 1)( −−−= 10 mm

( ) tFtF cos==

( ) gFzzktm

tz −+−= )(

( ) ( ) gtmtFzzktztm )()( −=+

( )temmtm 1.0

10 1)( −−−= 10 mm

( ) tFtF cos==

Ora simuliamo un sistema termico

Consideriamo un sistema termico rappresentato da un volume V circondato da unaparete e contenente un fluido

Sia Te(t) [K] la temperatura esterna alla parete, Tf(t) [K] la temperatura del fluido interno al volume, e q(t) [J/s] una sorgente di calore interna al volume. Sia Cf [J/K] la capacita termica del fluido, e sia Kie (J/K s) il coefficiente di scambio termico tra interno ed esterno.

( )tT f ( )tTe( )tq

( ) ( ) ( ) ( )( )tTtTKtqtTC feieff −+=

Per tradurre una equazione differenziale in termini di una combinazione tra blocchi Simulinksi deve esplicitare l’equazione differenziale rispetto alle derivate di ordine piu elevato.

( ) ( ) ( ) ( ) ( ) ( ) ( ) tTKtTKtqC

CtT fieeie

f −+=−+=11

Con un blocco integratore (Integrator) e con dei blocchi “Constant” e “Gain” si può realizzareil seguente modello di simulazione che ipotizza dei valori costanti per q e Te

Si noti come all’ingresso dell’integratore, punto cui corrisponde il segnale dTf/dt, venga

“costruita” elemento per elemento la formula ( ) ( ) ( ) ( )( ) tTtTKtqC

tT feie

f −+=1

Per eseguire la simulazione bisogna assegnare un valore ai parametri Cf [J/K] e Kie [J/Ks],

alla condizione iniziale Tf(0), ed al valore costante dei segnali Te e q.

La condizione iniziale dell’integratore si specifica nella relativa finestra di dialogo

Tf_zero=298.16; %[K], pari a 25°C

q=2000; % J/s q=2kW

Te=323.16; %[K], pari a 50°C;

Cf=2e3; % J/K per 1 kg di sostanza

Kie=1e2; % J / K s , per 1 m^2 di superficie di scambio

Script di parametrizzazione

Si desidera visualizzare Tf sia in gradi Kelvin che in gradi centigradi in 2 finestre Scope separate.La conversione da °K a °C può essere realizzata con il blocco Fcn (libreria User Defined Functions

Si desidera anche modificare il profilo di q(t)

Si può utilizzare il blocco Step, libreria Sources.

Il blocco Fcn consente di realizzare una funzione statica tra un parametro di ingresso (scalare o vettoriale) e un parametro di uscita scalare.Si deve scrivere l’espressione della funzione, denotando con u (parola riservata) la variabile in ingresso al blocco

Si realizzi lo schema seguente

Grafico della temperatura Tf in gradi centigradi

clear allclose allclc

if isempty(find_system('Name','sist_termico'))open_system('sist_termico')

Tf_zero=298.16; %[K], pari a 25°Cq=3000; % J/s q=2kWTstop=120; %secondi

Te=303.16; %[K], pari a 30°C; Cf=2e3; % J/K per 1 kg di sostanzaKie=1e2; % J / K s , per 1 m^2 di superficie di scambio

sim('sist_termico');figure(1)plot(tempo,Tf),gridTf_3000=Tf;

q=3500;sim('sist_termico');figure(2)plot(tempo,Tf),gridTf_3500=Tf;

q=4000;sim('sist_termico');figure(3)plot(tempo,Tf),gridTf_4000=Tf;

figure(4)plot(tempo,Tf_3000,tempo,Tf_3500,tempo,Tf_4000),gridlegend('q=3 kW','q=3.5kW','q=4kW')

Esecuzione automatica di test

Il seguente codice lancia in sequenza tre simulazioni con valori diversi di q

IMP. Utilizzo della funzione sim(‘model’)

Realizziamo lo stesso modello in maniera più compatta, utilizzando il blocco “Transfer Function” (libreria Continuous).

( ) ( ) ( ) ( )sTKsQsTKsTsC eiefieff +=+ ( ) ( ) ( )sTKsQsTKsC eiefief +=+

( ) ( ) ( )( )tTKtqKsC

sT eie

Notazione impropria, ma chiara

ief KsC +

( )tTK eie

Schema Simulink compattoSchema a blocchi

( ) ( ) ( ) ( )( ) ( ) ( ) ( )tTKtTKtqtTtTKtqtTC iieeieieieff −+=−+=

( ) ( ) ( )( )sTKsQKsC

sT eie

Trasformiamo secondo Laplace l’equazione differenziale

( )tT f

Parametrizzazione del blocco Transfer Function

Si devono specificare i coefficienti dei polinomi anumeratore e denominatore della FdT utilizzando lanotazione Matlab per la rappresentazione dei polinomi(un vettore che contiene i coefficienti dl polinomio inordine decrescente rispetto alle potenze di s)

( )ief KsC

Il blocco deve rappresentare la Funz. di Trasf.

] [ iefief KCKsC +

( ) ( ) ( )tutTKtTC fieff =+

( ) ( ) ( )tTKtqtu eie+= ] [ ief KC

]1[( )tu ( )tT f

Si può “aggirare” la trasformazione secondo Laplace ragionando sui coefficienti della equazione differenziali

][ 01 bb( )tu ( )tT

( ) ( ) ( ) ( ) ( )tubtubtTatTatTa 01012 +=++

012 aaa

FdT del secondo ordine

( ) ( ) ( ) ( ) ( )tubtubtyatyatya 01012 +=++

Far variare nel tempo la temperatura esterna ed osservare le corrispondenti fluttuazionidella temperatura interna

Nell’ipotesi che q(t) possa assumere solo due valori (il valore nullo q=0, ed un valorecostante q=Q* ) realizzare un sistema di controllo ON-OFF per regolare ad un valoredesiderato la temperatura del fluido interno al volume

Confrontare i risultati ottenuti impiegando come controllore un rele’ con e senza isteresi

Obbiettivi aggiuntivi

Schema SIMULINK FILE: termico1dof_feedback.mdl

Tf_zero=298;.16 %[K], pari a 25°C

q=2000; % J/s

Te=323;.16 %[K], pari a 50°C;

Cf=2e3; % J/K

Kie=1e2; % J / K s

Sistema a ciclo aperto

Ampiezzaisteresi = 4°C

Ampiezzaisteresi = 2°C

Complichiamo il modello

Trattiamo in maniera distinta gli accumuli termici nel volume e nella parete di contorno. Definiamo quindi un modello più complesso che mi fornisca anche l’evoluzione temporale della temperatura della parete Tp(t).

Cf [J/K] è la capacita termica del gas interno al volumeCp [J/K] è la capacita termica del materiale che costituisce la pareteKip [J/K s] è il coefficiente di scambio termico tra l’interno del volume e la parete. Kpe [J/K s] è il coefficiente di scambio termico tra la parete e l’esterno.

( )tT f ( )tTe( )tq

( )tTptemperatura della parete

( ) ( ) ( ) ( )( )tTtTKtqtTC fpipff −+=

( ) ( ) ( )( ) ( ) ( )( )tTtTKtTtTKtTC fpippepepp −−−=

Esplicitiamo il sistema di equazioni rispetto alle derivate di ordine più elevato

( ) ( ) ( ) ( )( )tTtTC

CtT fp

f −+=1

Importiamo nella pagina di lavoro un ulteriore blocco integratore la cui uscita sarà la temperatura Tp della parete.

Utilizziamo anche un visualizzatore a Display (libreria Sinks)Lo schema può essere realizzato come segue

( ) ( ) ( )( ) ( ) ( )( )tTtTC

KtT fp

p −−−=

Script di configurazione dei parametri

Tf_zero=25+273.15; %[K], pari a 25°C

Tp_zero=10+273.15; %[K], pari a 10°C

Cf=2e3; % J/K

Cp=10e3; % J/K

Kip=1e2; % J/K s

Kpe=1e2; % J/K s

q=2000; % J/s

Te=20+273.15 %[K], pari a 20°C;

Risultati della simulazione

Si desidera acquisire il segnale q(t) da un file esterno.

File “dati_problema.dat”

Il segnale q(t) sia disponibile nella forma di un file dati di tipo ASCII

La prima riga riporta i tempi e la seconda riporta il valore del segnale

Bisogna creare un mat-file (lo chiamiamo dati_mat.mat) che contenga tali informazioni.

%clear all

load dati_problema.dat

M=dati_problema;

tempi=M(1,:); %non utilizzato

segnale=M(2,:); %non utilizzato

save dati_mat M %si puo’ fare direttamente >> save dati_mat dati_problema

Tf_zero=298.16; %[K], pari a 25°C

Tp_zero=10+273.16; %[K], pari a 10°C

Cf=2e3; % J/K

Cp=10e3; % J/K

Kip=1e2; % J/K s

Kpe=1e2; % J/K s

q=2000; % J/s

Te=20+273.16; %[K], pari a 20°C;

Per importare il segnale si può utilizzare il blocco Simulink “From File”, dalla libreria Sources

Schema Simulink

In accordo con i dati acquisiti, la durata della simulazione deve essere posta pari a 2.

Creazione di sottosistemi

Si vuole rendere più compatta la rappresentazione del modello attraverso la definizione di un macroblocco come in figura

( )tT f

( )tTe

( )tTp

Simulink consente di definire dei macroblocchi (sottosistemi) che rappresentano una particolare interconnessione tra altri blocchi, alla quale si accede esplorando il contenuto del sottosistema.

Creazione di sottosistemi

Si dispongano i blocchi del modello come in figura, in modo che sia possibile tracciare un rettangolo nel quale entrano i segnali di input del sottosistema e dal quale escono i segnali di output.

Con il mouse si deve “tracciare” nella pagina di lavoro un rettangolo come quello nella slide precedente. Poi, dal menu Edit della pagina di lavoro Simulink, si deve selezionare il comando Create Subsystem

Ridisponendo i blocchi dello schema si puo realizzare la seguente configurazione

Cliccando sul blocco Subsystem si accede al suo contenuto.Si notino i blocchi “In1”,

“In2” , “Out1” ed “Out2”

Rinominando i blocchi “In1”, “In2”,

“Out1” ed “Out2” si può fare inmodo che nelle porte di input edi Output del Subsystemcompaia il nome della grandezzaassociata

La parametrizzazione del modello è attualmente effettuata per mezzo del file di scriptprecedentemente illustrato, che assegna una valore alle costanti simboliche usate nel modelloSimulink definendo opportune variabili omologhe nel workspace di Matlab. Tale m-file discript deve sempre accompagnare il file Simulink (che ha estensione .mdl) e deve esserelanciato prima di quest’ultimo

Può essere conveniente disporre di una maschera di parametrizzazione interna al modello Simulink.

In questo modo tutte le informazioni associate al modello sono contenute in un unico file (il file .mdl) e la riparametrizzazione del modello avviene in maniera più semplice, senza dovere ogni volta rilanciare il file script con estensione .m

Creazione di MASK

Bisogna portarsi con il mouse sul Subsystem“Camera di Combustione”, premere il tasto destro del mouse, e selezionare MaskSubsystem dal menu che compare

Dopo avere selezionato Mask Subsystem si deve nuovamente portarsi con il mouse sul Subsystem “Camera di Combustione”, premere il tasto destro del mouse, e selezionare stavolta Edit Mask dal menu che compare.

Finestra di configurazione della mask

Andare nel sottomenù:Parameters

Nel sottomenù Parameters completare le voci “Dialog Parameters” come in Figura, creando una riga per ognuno dei 6 parametri da settare .

Pulsante “Add”

Per ciascun parametro , oltre al nome della relativa variabile (Cf, Cp, …) si può riportare una frase descrittiva che comparirà nella maschera di configurazione, che riportiamo nella slide seguente.

Aggiunge una riga.

Maschera di parametrizzazione del Subsystem

Viene visualizzata facendo doppio click sul blocco Simulink del Sottosistema

E’ possibile inserire manualmente il valore dei parametri nelle opportune caselle di testo. I valori inseriti vengono memorizzati al salvataggio del file, e riproposti alla sua riapertura.

%clear all

M=dati_problema;

save dati_mat M

q=2000; % J/s

Te=20+273.16; %[K], pari a 20°C;

Script semplificato

Per accedere al contenuto del blocco “mascherato” premere il tasto destro del mouse e poi selezionare “Look Under Mask”

termico_2dof_matfile_mask01.mdl

Si inserisca una ulteriore mask che consenta l’impostazione manuale di q e Te

termico_2dof_matfile_mask02.mdl

%clear all

M=dati_problema;

save dati_mat M

Script

Sospensioni attive

Modello “quarter car” a due masse

Sospensioni attive

Modello “quarter car” a due masse

kgM b 400= kgM t 50=

mNk /102 5= mNkt /105.2 5= mNsc /102 3=

V(t) = volume di liquido contenuto nel serbatoio (in m3)

( ) ( )tawHtV =H(t) = livello del liquido nel serbatoio (in m)

Serbatoio a sezione costante rettangolare

portate volumetriche in ingresso ed in uscita dal serbatoio (in m3/sec)

( )tqi

( )tqu

( )tqi ( )tqu

( ) ( ) ( )tqtqtV ui −= ( ) ( ) ( ) ( )tqtqtHawtV ui −==

( ) ( ) ( )tqaw

11−=

Serbatoi idraulici

( ) HCtqu =

V(t) = volume di liquido contenuto nel serbatoio (in m3)

H(t) = livello del liquido nel serbatoio (in m)

Serbatoio a sezione rettangolare variabile

portate volumetriche in ingresso ed in uscita dal serbatoio (in m3/s)

( )tqi

( )tqu

( )tqi ( )tqu

( ) ( ) ( )tqtqtV ui −=

portate volumetriche in ingresso ed in uscita dal serbatoio (in m3/sec)( )tqi ( )tqu

Eq. di bilancio ( ) HCtqu =

Serbatoio a sezione rettangolare variabile( )tqi

( )tqu

( ) ( ) ( )tcwHtHH

bwtV += 2

r ( ) ( ) wH

( ) ( ) ( )tHcwtHH

( )( )

( ) ( )( )tqtq

bwtH ui −

( ) HCtqu =

( )tqi1

( ) ( )tqtq iu 21 =

( )tqu2

Serbatoi accoppiati in configurazioneverticale

( )( )

( ) ( )( )tqtq

wbtH ui 22

( ) ( ) ( )tqwa

tH ui 1

11−=

( ) ( ) 1121 HCtqtq iu ==

( ) 222 HCtqu =

( )tq1

( )tq3

( )tq12

( )tq2

( )tq23

( )( )

( ) ( ) ( )( )tqtqtq

wbtH 23122

( ) ( ) ( ) tqtqwa

tH 121

( ) ( )21211212 sgn HHHHCtq −−=

( ) ( )32322323 sgn HHHHCtq −−=

( ) ( ) ( ) tqtqwa

tH 233

Serbatoi accoppiati in configurazioneorizzontale

Setting Up Your Compiler

Before building your Embedded MATLAB Function block, you must set up your C compiler by running the mex -setup command, as described in the documentation for mex in the MATLAB Function Reference. You must run this command even if you use the default C compiler that comes with the MATLAB product for Microsoft Windows platforms. You can also use mex to choose and configure a different C compiler, as described in Building MEX-Files in the MATLAB External Interfaces documentation.

Supported Compilers for Simulation Builds. To view a list of compilers for building models containing Embedded MATLAB Function blocks for simulation:

Navigate to the Supported and Compatible Compilers Web page

[http://www.mathworks.com/support/compilers/current_release/]

Select your platform.

In the table for Simulink and related products, find the compilers checked in the column titled Simulink (Embedded MATLAB).

Esercizio

X0=0.1;Y0=0.1;

A=1;B=1;C=5;D=1;

xq ( )

BxyAxqM

( )qMq =

Modello in forma matriciale

function [ out ] = M_LV( in )

global A B C D

X=in(1); Y=in(2);

out(1)=A*X-B*X*Y;

out(2)=-C*Y+D*X*Y;

Function file “M_LV.m” da creare nella Current Directory

clear all

global A B C D

q_zero=[0.1 0.1];

( ) ( ) ( ) ( )txatxatxatx nn12121111 ...++=

( ) ( ) ( ) ( )txatxatxatx nn22221212 ...++=

( ) ( ) ( ) ( )txatxatxatx nnnnnn ...2211 ++=

( ) ( )tt Axx =

Sistema omogeneo di equazioni differenziali lineari

A=diag(-[1 2 3 4 5]);

B=rand(5);

A=inv(B)*A*B;

xzero=[1 1 1 1 1];

( ) ( ) ( ) ( ) ( ) ( ) ( )tubtubtubtxatxatxatx mmnn 121211112121111 ...... +++++++=

( ) ( ) ( ) ( ) ( ) ( ) ( )tubtubtubtxatxatxatx mmnn 222212122221212 ...... ++++++=

( ) ( ) ( ) ( ) ( ) ( ) ( )tubtubtubtxatxatxatx mnmnnnnnnnn +++++++= ...... 22112211

Sistema non omogeneo di equazioni differenziali lineari

A=diag(-[1 2 3 4 5]);

rand=rand(5);

A=inv(rand)*A*rand;

B=[1 0;0 1;1 1;0 0;1 2]

xzero=[1 1 1 1 1];

( ) ( ) ( )tBtt uAxx +=

2=m 5=n

( ) ( )txtx 21 =

( ) ( ) ( ) ( ) ( )ttxtxtxtx sin432 3212 +++=

( ) ( ) ( ) ( ) ( )tttxtxtx cossin4 313 ++−=

( )( )

( ) 101 =x

( ) 202 =x

( ) 303 =x

( ) ( ) ( )tBtt uAxx +=

( )( ) ( )

cossin

Rappresentazioniequivalenti

clear allclcA=[0 1 0;2 3 4;1 0 -1];B=[0 0;1 0;1 1];xiniziale=[1 2 3];

( )( )

Embedded Matlab Function block

“Matlab Function” block in Matlab R2013

Consente di eseguire un “function file” Matlab direttamenteall’interno di un modello Simulink

Default

Se si fa doppio click sul blocco lo si apre nell’editor

Modifichiamo il file

function [y1,y2] = fcn(u1, u2, u3)%#codegen

y1 = sin(u1)+u2;y2= u1*u2+u3

Tre ingressi e due uscite

L’aspetto del blocco cambia e compaiono i nuovi terminali diinput ed output

function [y1,y2] = fcn(u1, u2, u3)%#codegen

y1 = sin(u1)+u2;y2= u1*u2+u3

All’interno di un Matlab function non sono visibili le variabili delworkspace, neanche se queste vengono definite come globali

Embedded Matlab Function block

Consideriamo nuovamente il sistematermico del secondo ordine

( ) ( ) ( ) ( )( )tTtTC

CtT fp

f −+=1

( ) ( ) ( )( ) ( ) ( )( )tTtTC

KtT fp

p −−−=

( )tT f ( )tTe( )tq

( )tTptemperatura della parete

FILE: termico_2dof_matfile_mask01_EMF.mdl

Realizziamo il modello SIMULINK in maniera differente, generando i segnali q e Te con deiblocchi Signal Builder (libreria Sources) …..

….. e realizzando le equazioni per mezzo di un blocco Embeddbed MATLAB Function (EMF)

Integratore“vettoriale”

(bidimensionale)

y = usciteu = ingressip = parametri

Il blocco Camera di Combustione conserva la medesima maschera dellaimplementazione precedente

Codice del blocco EMF

function ydot = fcn(y,u,p)

Cf=p(1);

Cp=p(2);

Kip=p(3);

Kpe=p(4);

q=u(1);

Te=u(2);

Tf=y(1);

Tp=y(2);

Tfdot=q/Cf+(Kip/Cf)*(Tp-Tf);

Tpdot=(Kpe/Cp)*(Te-Tp)-(Kip/Cp)*(Tp-Tf);

ydot=[Tfdot Tpdot];Variabili globali non concesseall’interno di blocchi EMF.

I parametri della mask non sonodirettamente accessibili da partedella EMF

Codice analogo a quellodi un Function file

Parametrizzazione dei Signal Builder

Target

Zoom on T

Spostiamo verso sinistraquesta linea(doppio click su di essa)

Spostiamo verso destraquesta linea con modalitàanaloghe

Spostiamo verso l’alto questalinea con modalita analoghe

Profilo completo q(t) Profilo più complicato per q(t)

Proflo piu complicato per T(t)

Non serve piu un file script in abbinamento

Risultati della simulazione

Servomotore in corrente continua

( )( ) ( ) ( )tv

tdi V 1+−−=

( )( ) ( ) ( )tT

T 1−−=

Servomotore in corrente continua

( )( )( )

−−=

( )( ) ( ) ( )tv

tdi V 1+−−=

( )( ) ( ) ( )tT

T 1−−=

( ) ( ) ( )tBtt uAxx +=

SIMULAZIONE DI MANIPOLATORI ROBOTICI

Piano orizzontale

m1 , m2 masse dei link

l1 , l2 lunghezze dei link

J1 , J2 momenti di inerzia dei link

r1 , r2 coeff. di attrito viscoso

C1(t) , C2(t) coppie applicate ai giunti

Modello matematico

Termini inerziali

Attrito viscoso

Coppie applicate

Coriolis e centripete

1lmJlmm ++= 2

1lmJm +=

Coppie applicate

Siano le coordinate costanti di un punto di lavoro desideratodd

( ) ( ) ( )( )tKtKtC d

pd −+−= 1

Si applichino le seguenti coppie ai giunti (controllore PD)

( ) ( ) ( )( )tKtKtC d

pd −+−= 2

pd KK , “guadagni” costanti

( )( )

( )( )( )

( ) ( ) ( ) ( )( ) ( ) ( ) tqRtCtqtqFqMtq inNL −+= − ,1

Vettore delle variabili di giunto

Modello vettoriale in forma esplicita

m1star=0.25*m1*l1^2+J1+m2*l1^2;

m2star=J2+0.25*m2*l2^2;

q0=[pi/2;0];

q0dot=[0;0];

R=[r1 0;0 r2];

Kp=100;

Kd=20;

( ) ( ) ( ) ( )( ) ( ) ( ) tqRtCtqtqFqMtq inNL −+= − ,1

function [ out ] = generaA( in )

m1star=0.25*m1*l1^2+J1+m2*l2^2;

m2star=J2+0.25*m2*l2^2;

alfa=in(1);

beta=in(2);

out(1,1)= m1star+m2star+m2*l1*l2*cos(beta);

out(1,2)= m2star+0.5*m2*l1*l2*cos(beta) ;

out(2,1)=m2star+0.5*m2*l1*l2*cos(beta);

out(2,2)=m2star;

function [ out ] = generaFinNL( in )

m1star=0.25*m1*l1^2+J1+m2*l2^2;

m2star=J2+0.25*m2*l2^2;

alfa=in(1);

beta=in(2);

alfadot=in(3);

betadot=in(4);

out(1,1)=m2*l1*l2*(alfadot*betadot+0.5*betadot^2)*sin(beta);

out(2,1)=-m2*l1*l2*0.5*alfadot^2*sin(beta);

Function files

Modello compatto con Embedded Matlab Function

function qddot = fcn(q,qdot,C)

m1star=0.25*m1*l1^2+J1+m2*l2^2;

m2star=J2+0.25*m2*l2^2;

alfa=q(1);

beta=q(2);

alfadot=qdot(1);

betadot=qdot(2);

A=[m1star+m2star+m2*l1*l2*cos(beta)

m2star+0.5*m2*l1*l2*cos(beta);

m2star+0.5*m2*l1*l2*cos(beta) m2star];

FinNL=[m2*l1*l2*(alfadot*betadot+0.5*betadot^2)*sin(

beta);

-m2*l1*l2*0.5*alfadot^2*sin(beta)];

R=[r1 0;0 r2];

qddot = inv(A)*(FinNL+C-R*qdot);

Embedded Matlab Function

Manipolatore planare a 2 gdl

( ) ( ) ( ) ( ) ( ) =++ qgtqqBtqqH

( )qH ( )qB ( )qg

( ) ( ) ( ) ( ) ( ) =++ qgtqqBtqqH

( ) ( ) ( ) ( ) ( ) qgtqqBqHtq −−= − 1

( ) ( ) ( )qgtqqBF +=

function [ out ] = generaA_manip( in )

lc1=1;

lc2=1;

q1=in(1);

q2=in(2);

out(1,1)= m1*lc1^2+I1+m2*(l1^2+lc2^2+2*l1*lc2*cos(q2))+I2;

out(1,2)= m2*l1*lc2*cos(q2)+m2*lc2^2+I2;

out(2,1)=out(1,2);

out(2,2)=m2*lc2^2+I2;

function [ out ] = generag_Bqdot( in )

lc1=1;

lc2=1;

g=9.81;

q1=in(1);

q2=in(2);

q1dot=in(3);

q2dot=in(4);

qdot=[q1dot;q2dot];

h=m2*l1*lc2*sin(q2)

B=[-h*q2dot -h*(q1dot+q2dot); h*q1dot 0];

grav=[m1*lc1*g*cos(q1)+m2*g*(lc2*cos(q1+q2)+l1*cos(q1)); m2*lc2*g*cos(q1+q2)];

out=B*qdot+grav;

Sistemi termici spazialmente distribuiti

Variabile di interesse: la temperatura

Barra metallica filiforme di lunghezza L = 50 cm

0=x Lx =

( )txT ,

Ipotesi : barra uniforme, termicamente isolata, in assenza di sorgenti di calore

( ) ( )( )txTk

txTCp ,

Equazione del calore monodimensionale

= è il coefficiente di diffusione [m2/s]

pC è la capacita termica a pressione costante per unita di massa [J/kg K]

è la conduttività termica [J/s m K ; W / mK ]k

è la densità [kg /m3 ]

( ) ( )( )txTk

txTCp ,

( ) ( )2

Valori per il ferro puro (rif. Schaum Trasmissione del calore, Tabella B-1 pag. 306)

J/kgK 435=pc C) 100 (a W/mK 63 =kkg/m 7856=

C) 100 (a /m 108.1 25 = − s

Analisi in regime stazionario

0=x Lx =

h h h h h 5/Lh =

4 “NODI SOLUZIONE”

1xx = 2xx = 3xx =4xx =

La soluzione è calcolabile in forma chiusa, ma scegliamo di risolvere l’equazionenumericamente perché questo modo di ragionare di consentirà di risolvere agevolmenteanche i problemi distribuiti in regime transitorio , che vedremo in seguito.

( ) ( )xTtxT ,

Siano noti e mantenuti costanti i valori di temperatura agli estremi ( )00 TT = ( )LTTL =

( ) xL

TTTxT L

−+= 0

0Soluzione in forma chiusa

Definiamo il vettore che contiene le temperature incognite nei nodi soluzione.

TTTTT 4,321 .,,=T

Approssimiamo l’equazione in un intorno dei nodi soluzione, e mettiamo a sistema

02 012 =+− TTT

02 123 =+− TTT

02 234 =+− TTT

02 34 =+− TTTL

Isoliamo alla sinistra dell’uguale le temperature incognite, ordinandole

0212 TTT −=+−

02 321 =+− TTT

02 432 =+− TTT

LTTT −=− 43 2

Sistema di equazioni lineari

T0=10; %temperatura nel punto x=0TL=50; %temperatura nel punto x=LL=0.5; %Lunghezza della sbarraN=4; %numero di nodi di discretizzazioneh=L/(N+1);

A=[-2 1 0 0;1 -2 1 0; 0 1 -2 1; 0 0 1 -2];

b=[-T0;0;0;-TL];

disp('il vettore incognito è:')T=inv(A)*b

disp('Vettore completo:')T_tot=[T0 T' TL]'

x=0:h:L;T_teorica=T0+(TL-T0)/L*x

plot(x,T_tot,'*',x,T_teorica,'r'),gridxlabel('Ascissa lungo la sbarra [m]')ylabel('Temperatura [°C]')

0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.510

Ascissa lungo la sbarra [m]

Modificare il codice per includere un numeroarbitrario N di nodi soluzione.

T0=10; %temperatura nel punto x=0TL=50; %temperatura nel punto x=LL=0.5; %Lunghezza della sbarraN=10; %numero di nodi di discretizzazioneh=L/(N+1);

A=zeros(N,N);for i=1:N,

A(i,i)=-2;endfor j=1:N-1,

A(j,j+1)=1;A(j+1,j)=1;

b=zeros(N,1);b(1)=-T0;b(N)=-TL;

disp('il vettore incognito è:')T=inv(A)*b%ISTRUZIONE ALTERNATIVAT=A\b;

disp('Vettore completo:')T_tot=[T0 T' TL]'

x=0:h:L;T_teorica=T0+(TL-T0)/L*x

plot(x,T_tot,'*',x,T_teorica,'r'),gridxlabel('Ascissa lungo la sbarra [m]')ylabel('Temperatura [°C]')

0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.510

Ascissa lungo la sbarra [m]

La soluzione in regime permanente è ancora calcolabile in forma chiusa

E’ noto e prefissato il valore di temperatura T0 all’estremo sinistro, mentre latemperatura all’estremo destro non è imposta.

Ora studiamo un problema differente

Poiché la sbarra è termicamente isolata si avrà che il flusso termico nei bordodestro è pari a zero

( ) 0=

( ) 0TxT =Soluzione in forma chiusa

( ) ( )

−−

0212 TTT −=+−

02 321 =+− TTT

02 432 =+− TTT

043 =−TT

Modificare il codice Matlab

021 =+− TT

02 321 =+− TTT

02 432 =+− TTT

043 =−TT

Ora analizziamo il caso in cui anche la temperatura all’estremo sinistro non èimposta.

Il sistema lineare non puo essere risolto numericamente percheammette infinite soluzioni

( ) .cost=xTSoluzione in forma chiusa

Seguendo il medesimo ragionamento adottato in precedenza si ottiene il sistema:

Cosa succede se le due temperature agli estremi, ancorché impostedall’esterno, siano variabili nel tempo ?.

L’analisi in regime stazionario perde di significato !

( ) ( )2

L’equazione del calore deve ora essere considerata nella sua interezza.

Il ragionamento svolto in precedenza va opportunamente modificato.

0=x Lx =

h h h h h 5/Lh =

1xx = 2xx = 3xx =4xx =

La discretizzazione del dominio è sempre valida.

( ) ( ) ( ) ( ) ( ) TtTtTtTtTt 4321 ,,,=T

( )tT0 ( )tTL

Definiamo il vettore che contiene le temperature incognite nei nodi soluzione, che adifferenza dal caso stazionario non sono più delle costanti, ma sono delle funzioni del tempo.

Condizioni al contorno NOTE

( ) ( )2

( ) ( )tTtT 0,0 =

( ) ( )tTtLT L=,Questa tipologia di condizioni al contorno sono dette di “DIRICHLET”.

( ) ( )txTtT ii ,= 4,3,2,1=i5

Lihixi ==

( ) ( )

( )( ) ( ) ( )

tTtTtTtxT

( )( ) ( ) ( )

tTtTtTtxT

( )( ) ( ) ( )

tTtTtTtxT

( )( ) ( ) ( )

tTtTtTtxT

( )( )

2 , 4,3,2,1=i

Le derivate seconde vengono approssimate mediante differenze finite, in maniera analoga aquanto fatto in precedenza ma coinvolgendo funzioni del tempo anzichè costanti

( )( )

2 , 4,3,2,1=i

( ) ( ) ( ) ( )( )tTtTtTh

tT 01221 2 +−=

( ) ( ) ( ) ( )( )tTtTtTh

tT 12322 2 +−=

( ) ( ) ( ) ( )( )tTtTtTh

tT 23423 2 +−=

( ) ( ) ( ) ( )( )tTtTtTh

tT L 3424 2 +−=

Sostituendo le approssimazioni alle differenze finite si ottiene un sistema di equazionidifferenziali ordinarie (ODE), di ordine pari al numero di nodi soluzioni (4 nel caso in esame).

( ) ( ) ( )( ) ( )tTh

tT 022121 2

++−=

( ) ( ) ( ) ( )( )tTtTtTh

tT 32122 2 +−=

( ) ( ) ( ) ( )( )tTtTtTh

tT 43223 2 +−=

( ) ( ) ( )( ) ( )tTh

tT L24324 2

Risulta conveniente separare alla destra dell’uguale le quantità incognite dalle quantità noteT0(t) e TL(t)

( ) ( ) ( ) TN tTtTt 10 , +=u

( ) ( ) ( )ttt BuATT +=

( ) ( ) ( ) ( ) ( ) TtTtTtTtTt 4321 ,,,=T

Devono essere note le temperature nei nodi soluzione all’istante iniziale t=0, cioè il vettorecostante

Ora sia il coefficiente di diffusione che lo step h di discretizzazione spaziale influenzanoesplicitamente la soluzione. L’analisi in regime stazionario prevedeva la scomparsa di talicoefficienti dalle relazioni risolutive

( ) ( ) ( ) ( ) ( ) TTTTT 0,0,0,00 4321=T

Va rimarcato che trascorso un transitorio sufficientemente lungo i profili di temperatura nondipendono più dalle condizioni iniziali, e tendono ad una soluzione di “regime dinamico” lecui caratteristiche sono indipendenti dalle condizioni iniziali

clear allclc

L=0.5;N=20;

rho=7856;Cp=435;k=63;alfa=k/(rho*Cp);

h=L/(N+1);

A(j,j+1)=1;A(j+1,j)=1;

B=zeros(N,2);B(1,1)=1;B(N,2)=1;

Ac=(alfa/h^2)*A;Bc=(alfa/h^2)*B;

Tcamp=0.1;Tfinale=100;

x=linspace(h,L-h,N);t=0:Tcamp:Tfinale;numsample=length(t);

Tinit=10*ones(1,N);

T_soluz=zeros(numsample,N);T_dot=zeros(numsample,1);

T_soluz(1,:)=Tinit;

T0=30+5*sin(0.1*t)';TL=30*ones(numsample,1);

for i=1:numsample-1Tdot=Ac*T_soluz(i,:)'+Bc*[T0(i);TL(i)];T_soluz(i+1,:)=T_soluz(i,:)+Tcamp*Tdot';

x_ext=linspace(0,L,N+2);T_soluz_ext=[T0 T_soluz TL];

[X,Y] = meshgrid(x_ext,t);

h=mesh(X,Y,T_soluz_ext);title('Distribuzione temperatura lungo la sbarra.');xlabel('Coordinata spaziale x [m] ','FontName','times','FontSize',14);ylabel('Tempo [s]','FontName','times','FontSize',14);zlabel('T(x,t)','FontName','times','FontSize',14);set(gca,'FontSize',14,'FontName','Times');

Eulero esplicito

Cond. al contorno

Grafica

Profilo spaziotemporale della soluzione

clear allclc

L=0.5;N=5;

rho=7856;Cp=435;k=63;

alfa=k/(rho*Cp);

h=L/(N+1);

A(j,j+1)=1;A(j+1,j)=1;

B=zeros(N,2);B(1,1)=1;B(N,2)=1;

Tinit=30*ones(1,N);

T_soluz(1,:)=Tinit;

x_ext=linspace(0,L,N+2);T_soluz_ext=[T0 T_soluz TL];

( ) ( )2

( ) ( )tTtT 0,0 =

( ) 0, =

TQuesta tipologia di condizione al contorno viene detta di “NEUMANN”.

Ora applichiamo una condizione al contorno che non “fissi” il valore della temperaturanell’estremo destro, ma imponga l’isolamento termico dell’estremo destro x=L inanalogia con quanto fatto nella analisi a regime.

Si ricava il modello discretizzato( )( ) ( )

tTtTtL

T NL −

( ) ( ) TtTt 0 , 0=u

( ) ( ) ( )ttt BuATT +=

( ) ( ) ( ) ( ) ( ) TtTtTtTtTt 4321 ,,,=T

( ) ( )tTt 0=u

( ) ( ) ( )ttt BuATT +=

( ) ( ) ( ) ( ) ( ) TtTtTtTtTt 4321 ,,,=T

clear allclc

L=0.5;N=20;

rho=7856;Cp=435;k=63;alfa=k/(rho*Cp);

h=L/(N+1);

A=zeros(N,N);

for i=1:N,A(i,i)=-2;

endfor j=1:N-1,

A(j,j+1)=1;A(j+1,j)=1;

endA(N,N)=-1;

B=zeros(N,2);B(1,1)=1;B(N,2)=1;

Tinit=10*ones(1,N);

T_soluz(1,:)=Tinit;

x_ext=linspace(0,L,N+2);T_soluz_ext=[T0 T_soluz T_soluz(:,end)];

Sistema LTI (Linear Time Invariant)MIMO (multi-input-multi-output)

( )tT( )tT

( )tAT

( )tBu( )tu

( ) ( ) ( )ttt BuATT +=

Uso di Matrix Gain Integratore saturato

Un processo termico 2D a parametri distribuiti

Cilindro cavo.

Parti grigie in acciaio.

Parte bianca: volume con vapore ad alta temperatura.

Parte nera: piccolo volume interno

Es. Sezione di una turbina a vapore

statore

rotorerCoordinata radiale

Es. Sezione di una turbina a vapore

Ipotesi: simmetria angolare della distribuzione di temperatura

( )trT ,

Si desidera calcolare la distribuzione di temperatura nel rotore (in uno dei suoi raggi)

maxmin rrr

Mediante misure acquisite in una turbina in esercizio, si suppone nota la temperatura nella “parte bianca” (regione del vapore in alta temperatura)

Problema complesso (anche nella formulazione semplificata sotto esame) perche il modello matematico è una equazione alle derivate parziali (sistema a parametri distribuiti, sistema infinito-dimensionale)

In una modellazione più dettagliata, si potrebbe essere interessati a valutare la distribuzione di temperatura nella superficie delle pale di rotore, onde valutare gli stress termici sui materiali.

Eq. di diffusione (Equazione del calore) monodimensionale in coord. cilindriche , con unica variabile spaziale la coordinata cilindrica radiale r

= è il coefficiente di diffusione [m2/s]

pC è la capacita termica a pressione costante per unita di massa [J/kg K]

è la conduttività termica [J/K s m]k

è la densità [g /m3 ]

minrr =

maxrr =

hirri += min

1..., ,2 ,1 ,0 += Ni

minmax

( ) ( )trTtT ii ,=

Approssimazione delle derivate spaziali mediante differenze finite

( )( ) ( )

tTtTtT

( )( ) ( ) ( )

tTtTtTtT

−+ +−=

Sistema di ODE

( )( ) ( ) ( ) ( ) ( )

( ) ( ) ( )

( ) ( ) ( ) =

tTtTtT

Ni ..., ,2 ,1 =

Sistema di ODE

( ) ( ) ( ) ( )

+= tTtT

KtT 01

( ) ( ) ( ) ( )

+= tTtT

KtT 12

( ) ( ) ( ) ( )

+= −+ tTtT

KtT NN

N 11221

( ) ( ) ( ) TN tTtTt 10 , +=u

( ) ( ) ( ) ( ) TNN tTTTtTtTt ,...,,, 2,321 −=T

( )tT0 ( )tTN 1+ Boundary conditions ( ) ( ) ( )ttt BuATT +=

210000

000121

000012

( ) ( ) ( )ttt BuATT +=

( ) ( ) ( ) ( ) TNN TTTTT 0,...,,0,00 12,321 −−=T Condizioni iniziali

Modello Simulink

Due diverse modalità di specificare i due segnali di input.

Linee spesse che rappresentano segnali multidimensionali

Contenuto del Subsystem “Modello termico Rotore”

TurbinaVapore01.mdl

T_centrorotore.mat

T_vapor.mat

Esportazione dati verso il Workspace

Parametrizzazione del modello

Utilizziamo una mask

Costruzione della maschera di parametrizzazione

Calcolo delle Matrici A e B

Matrici A e B (per n = 6)

r_min=0.05;

r_max=0.392;

rho=7900;

C=0.45;

T_in=80;

h=(r_max-r_min)/n;

K=k/(rho*C*1e3)/h^2;

A=zeros(n,n);B=zeros(n,2);

for i=1:n,

A(i,i)=-(2+h/(r_min+i*h));

for j=1:n-1,

A(j,j+1)=1+h/(r_min+j*h);

A(j+1,j)=1;

B(1,1)=1;B(n,2)=1+h/r_max;

A=K*A; B=K*B;

Codice copiabileed seguibile

Documentation - descrizione della MASK

Profili temporali degli elementi del vettore T Uscita del blocco Scope “T”

Analisi dei risultati

Impostiamo preliminarmente due valori costanti per le temperature al contorno

( ) ( ) ( ) CtTtTtT

Profilo di temperatura al nodo 40

Con dei profili differenti per le boundary conditions, l’evoluzione del profilo di temperatura èdifferente.

Ora processiamo in Matlab i risultati della simulazione, creando dei grafici 3D.

r_min=0.05;

r_max=0.392;

[X,Y] = meshgrid(linspace(r_min,r_max,n),tout);

h=mesh(X,Y,Trot)

title('Distribuzione temperatura rotore.')

xlabel('Coordinata radiale r [m]

','FontName','times','FontSize',14)

ylabel('Tempo [s]','FontName','times','FontSize',14)

zlabel('T(r,t)','FontName','times','FontSize',14)

set(gca,'FontSize',14,'FontName','Times')

Invertiamo la direzione dell’asse dei tempi

Completiamo questo esempio mostrando l’impiego, nella Mask, di variabili di configurazionetipo “popup” o “checkbox” e un loro possibile impiego.

Si apportino le seguenti modifiche alla lista dei Parameters

Si apportino le seguenti modifiche alle istruzioni di Initialization

Si apportino le seguenti modifiche allo schema Simulink

Ora si esegua il modello con diverse scelte per le variabili pop up e checkbox, e si analizzino irisultati.

FILES:

TurbinaVapore02.mdl

T_centrorotore.mat

T_vapor.mat

SimMechanics

Modellazione di sistemi meccanici multi-body

Toolbox avanzati

Modello VMRL 3D

Vista differente

Variare il controllo sul FUEL (Subsystem Control)

Osservare la velocità dell’albero motore in RPM, e correlarne le variazioni ai corrispondenti incrementi e decrementi del fuel rate

Toolbox avanzati

SimDriveline

Componenti e modelliarea automotive

Animazione 3D

MATLAB-SIMULINK€¦ · 3 Finestra di avvio (v. 7.8.0) Cartella corrente Avvio SIMULINK Editor...

Documents

Transcript of MATLAB-SIMULINK€¦ · 3 Finestra di avvio (v. 7.8.0) Cartella corrente Avvio SIMULINK Editor...

BigMat Pisano: Flyer Febbraio 2014

Giugno Pisano 2016

Laboratorio di Matlabcontrol.dii.unisi.it/~control/MatLab/febbraio2009/...DII - Universtita di Siena` 23 Simulink Simulink e un pacchetto software di Matlab che permette di simulare

Guida Introduttiva a SIMULINK · 2019-10-04 · Guida introduttiva a SIMULINK dic a 7 Esempio introduttivo: visualizzazione di una sinusoide Come primo passo si apra una nuova pagina

Guida Introduttiva a SIMULINK

Studio in ambiente Matlab/Simulink del sistema di controllo di ...

01.03 Introduzione a Simulink

STUDIO IN MATLAB/SIMULINK DI UN SISTEMA DI ...tesi.cab.unipd.it/44160/1/TESI-[Bergamo_Elena_614521].pdfsimulazione del sistema attraverso Simulink. 4 5 Ringrazio il mio relatore, gli

Esercizi in MATLAB-SIMULINK Matlab.pdf · trovarlo nel Simulink Library Browser seguendo il percorso Simulink Extras\ Additional Linear\PID Controller. Affinché le specifiche del

Eclissi catalogo edoardo pisano

Biagio Pisano

Ing. Alessandro Pisano - Home - people.unica.it · 7 Appunti di Controlli Automatici – Sistemi con ritardo – v. 1.1 Ing. Alessandro Pisano – pisano@diee.unica.it Figura 5 –

INTRODUZIONE A SIMULINK - UnivAQing.univaq.it/digennaro/Didattica/CCAA/Materiale/Simulink.Manuale.pdf · Introduzione a SIMULINK 0 25 Lab Trascinare 2 blocchi Gain Math Per ruotare

Intro Simulink

Introduzione a Matlab/Simulink

Simulink - UNISAcesposito/materiale/lezioni/Argomento15_Simulink.pdf · Per utilizzare Simulink bisogna lanciarne l’interfaccia, digitandosimulink.Oppurecliccarel’appositaicona.

MATLAB-SIMULINK Ing. Alessandro Pisano pisano@diee.unica.it 1 Bergamo, 17-19/11/2009 Matlab.

Polo Ambientale Monte Pisano eventi escursioni attività 2013

E1. IMPLEMENTAZIONE in MATLAB-SIMULINK del ......1 E1. IMPLEMENTAZIONE in MATLAB-SIMULINK del MODELLO e del CONTROLLO di un MOTORE ELETTRICO a CORRENTE CONTINUA 1. EQUAZIONI DEL MODELLO

SISTEMI DI CONTROLLO Ingegneria Meccanica e Ingegneria …lbiagiotti/MaterialeDidattico0809/CA-02M... · Luigi Biagiotti Sistemi di Controllo Simulink intro -- 5 Simulink • Nelle