“Elevação e Estabilização Vertical do Pêndulo Invertido ...dulley/dulley_docs/Pendulo_Invertido_Rotacional... · 2.1 Modelagem Dinâmica ... Além da análise de modelagem

Escola Politécnica da Universidade de São Paulo

PTC 2531 Laboratório de Projeto de Automação II

“Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional”

(Relatório Final)

Orientadores: José Jaime da Cruz e Gilberto Martha de Souza

Integrantes da Equipe: NOME NºUSP Lucas Peetz Dulley 2933561 Germano Claro Simões Machado 2997570

Dezembro de 2003

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional

Índice 1 Introdução ................................................................................................................... 3

1.1 Objetivo............................................................................................................... 3 1.2 Visão Geral do Projeto........................................................................................ 3 1.3 Agradecimentos .................................................................................................. 4

2 O Pêndulo Invertido Rotacional.................................................................................. 5 2.1 Modelagem Dinâmica......................................................................................... 5 2.2 Linearização em torno do equilíbrio instável ( πθ =2 ).................................... 13

2.2.1 Estratégia 1 : Linearização entrada-saída ................................................. 13 2.2.2 Estratégia 2 : Linearização por Taylor...................................................... 17

2.3 Linearização em torno do equilíbrio estável ( 02 =θ ) ...................................... 22 3 Componentes Utilizados ........................................................................................... 25

3.1 Elétricos ............................................................................................................ 25 3.2 Mecânicos ......................................................................................................... 25

4 Identificação de Parâmetros...................................................................................... 29 4.1 Parâmetros do Motor......................................................................................... 29 4.2 Parâmetros do Sistema...................................................................................... 34

4.2.1 Momentos de Inércia................................................................................. 34 4.2.2 Coeficientes de Atrito Viscoso ................................................................. 36 4.2.3 Coeficientes de Atrito Estático ................................................................. 41

5 Controle..................................................................................................................... 43 5.1 Elevação............................................................................................................ 43 5.2 Equilíbrio .......................................................................................................... 46

5.2.1 Controle Linear em Tempo Contínuo ....................................................... 46 5.2.2 Controle Linear em Tempo Discreto ........................................................ 51

6 Implementação e Testes............................................................................................ 53 7 Conclusões ................................................................................................................ 58 8 Cronogama................................................................................................................ 59 9 Bibliografia ............................................................................................................... 60 10 Anexo I (Scripts de Matlab e Códigos Fonte) ...................................................... 61

10.1 Modelo Não Linear do Pêndulo Invertido (sistemaEuler.m)............................ 61 10.2 Modelo Linear: Funções de Transferência e Representação em Espaço de Estados (snormal_V.m)................................................................................................. 62 10.3 Elevação do Pêndulo (sistemaEuler.m) ............................................................ 63 10.4 Código Fonte do Controlador em C (lqrdig.c).................................................. 64

11 Anexo II (Desenho das Peças) .............................................................................. 68

2


1 Introdução 1.1 Objetivo O propósito deste relatório é exibir os procedimentos necessários à construção de um pêndulo invertido rotacional e, após construído, como elevá-lo e equilibrá-lo na posição vertical. Este documento contém uma análise teórica de modelagem de sistemas compostos por ligamentos, seguida de diversas abordagens escolhidas para os modos de controle do pêndulo. Além da análise de modelagem e controle, estão presentes todas as etapas da construção e identificação dos parâmetros do sistema. Os resultados simulados e experimentais estão comparados no final do relatório. 1.2 Visão Geral do Projeto Os dois pontos principais de todo o trabalho são a elevação e estabilização do pêndulo invertido. Naturalmente, tão importante quanto os algoritmos de controle empregados é a implementação destes em ambiente Matlab (Simulink). A rotina de elevação do pêndulo é baseada na estratégia de verificação de velocidades e posições de forma a sempre aumentar a energia do sistema. Uma vez elevado o pêndulo, o controle de equilíbrio passa a atuar, garantindo que pequenas perturbações não causem a queda do mesmo. Para esta etapa, foram testadas diversas técnicas de controle, desde controle linear por realimentação de estados (LQR –Linear Quadratic Regulator), passando por Modos Deslizantes até controle robusto por H infinito. A linearização do sistema foi feita empregando a técnica de Taylor, porém o cancelamento das linearidades utilizando realimentação entrada-saída também será discutido neste trabalho. Na etapa da construção mecânica serão indicados as principais considerações quanto à funcionalidade, facilidade de montagem-desmontagem, problemas enfrentados e sugestões para futuras implementações por parte de alunos de engenharia elétrica. Montado o sistema, seguiu-se a identificação e cálculo dos parâmetros para, enfim, serem aplicados nas simulações e nas execuções em laboratório.

3

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional 1.3 Agradecimentos Gostaríamos de agradecer profundamente às seguintes pessoas pela gentileza de nos ajudar em diversas etapas do projeto. Todas foram de extrema importância para o bom andamento dos trabalhos. - Professor José Jaime da Cruz : pela ajuda providencial na resolução de problemas que pareciam insolúveis, e por sempre se mostrar solícito e prestativo. - Professor Gilberto Martha de Souza: certamente, sem sua ajuda, a execução do nosso projeto teria se limitado às simulações. Sua contribuição para o desenvolvimento dos componentes mecânicos fez toda a diferença. - Amigos da oficina da engenharia mecânica: Cícero, Fidel, Laércio, Paulo e Edson. A paciência de vocês com nossos conhecimentos limitados de projeto mecânico e o socorro que nos prestaram na construção de cada componente os tornaram elementos fundamentais para que tudo desse certo. Obrigado mais uma vez pelos meses de ajuda.

4


2 O Pêndulo Invertido Rotacional 2.1 Modelagem Dinâmica O método de modelagem escolhido foi o de Newton-Euler, por ser iterativo e de cálculos mais simples (manipulações algébricas) do que o método de Lagrange (comumente utilizado para problemas desta natureza). Os ligamentos do manipulador são descritos através dos Parâmetros de Denavit-Hartenberg (D.H.). Também são considerados como corpos rígidos, desprezando deformações. Foram escolhidos três sistemas de coordenadas fixados aos ligamentos, sendo que o sistema 0 base é coincidente com a posição de 1, quando 01 =θ . A partir dos parâmetros de D.H. as matrizes de transformação dos sistemas de coordenadas foram calculadas.

A seguir ilustramos o desenvolvimento matemático passo a passo.

Vista Lateral

1z

1x1y 2z

2x

2y

Vista Frontal

1L

21 ˆˆ zy =2θ

1x

2x

2y

0x 1θ

10 ˆˆ zz =

1x

0y

1y

Vista Superior

Parâmetros Denavit-Hartenberg ai-1 dii αi-1 θi-1

1 0 0 0 θ1

2 -90º 0 L1 θ2

5

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional Parâmetros físicos do manipulador:

m1: massa do ligamento 1

nto 1

dotando as coordenadas escolhidas, obtemos a seguinte tabela, pela Teoria de

αi-1 ai-1 di θi-1

m2: massa do ligamento 2 L1: comprimento do ligameL2: comprimento do ligamento 2

AManipuladores: I 1 0 0 0 θ2 -90º 0 L1 θ2

Tabela 1 : parâmetros dos ligamentos

s matrizes de rotação-translação do sistema 2 em relação ao 1 e do sistema 1 em relação

As velocidades angular e linear do ligamento i+1 relativas ao ligamento i e expressas no

- velocidade angular :

velocidade linear:

Aao 0 são:

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡ −

=

1000010000cos00cos

11

11

01

θθθθ

sensen

T⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡ −=

1000cos0cos

11

1101 θθ

θθsen

senR

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

−−

−

=

100000cos

10000cos

22

1

22

12 θθ

θθ

senL

sen

T⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−−

−=

0cos1000cos

22

2212

θθ

θθ

sen

senR

sistema i+1 são :

1

11

11

1 ˆ.. ++

++

++ += i

iii

iiii

i zwRw θ&

-

6


).( 11

11

++

++ Ο∧+= i

ii

ii

iiii

i wvRv

abemos ainda que :

rocedemos com o cálculo das velocidades angulares

S

Tii

ii RR )( 1

1+

+ =

: P

[ ]TzwRw 11

110

0101

1 00ˆ.. θθ && =+=

(velocidade angular do centro de massa da barra vertical)

a seqüência, obtemos as velocidades lineares

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−

−

=⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡+

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−−−

=+=

2

21

21

2

1

22

22

22

2112

122 cos.

.

100

.00

.010

cos00cos

ˆ..

θ

θθ

θθ

θ

θ

θθθθ

θ&

&

&

&

&

&

sen

sensen

zwRw

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−

−

==

2

21

21

22

22 cos.

.

θ

θθ

θθ

&

&

& sen

wwC

: N

[ ]TwvRv 000).( 1

00

00

0101

1 =Ο∧+=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−

=Ο∧+=0

.

cos..

).( 211

211

21

11

112

122 θθ

θθ

senL

L

wvRv &

&

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

−

−−

=

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

∧⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−

−

+⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−

=∧+=

22

1

211

22211

2

2

21

21

211

211

22

22

22

22

2.

.2

.cos..

02

0

cos.

.

0.

cos..

θθ

θθ

θθθ

θ

θθ

θθ

θθ

θθ

senLsenL

LLL

sen

senL

L

CwvvC

&

&

&&

&

&

&

&

&

elo método de Newton-Euler, as demais grandezas (acelerações, forças, torques etc.) nas juntas são determinadas por equações iterativas:

P

7

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional -aceleração angular :

11

111

11

11 ˆ.)ˆ..( +

+++

++

++

+ +∧+= ii

iii

iii

iii

iii zzwwRw θθ &&&&&

celeração linear total ( radial + tangencial ) :

aceleração linear do centro de massa do ligamento i+1 :

ara levarmos em conta o efeito da gravidade, consideramos que a base O se ovimenta para cima com uma aceleração g.

-a

])().[( 111

11

ii

ii

ii

ii

ii

iii

iii vwwwRv &&& +Ο∧∧+Ο∧= ++

++

+

-

11

11

11

11

11

11

11 )( +

++

++

++

++

++

++

+ +Ο∧∧+Ο∧= ii

Ciii

ii

ii

Ciii

ii

Ciii vwwwv &&&

PmOs resultados calculados foram : -acelerações angulares:

[ ]Tw 111 00 θ&&& =

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−

−

==

2

21

21

22

22 cos.

.

θ

θθ

θθ

&&

&&

&&

&&

sen

wwC

acelerações lineares- :

, onde

[ ]Tgv −= 0000 & [ ]Tg 78.900 −=

[ ]Tgv −= 0011 &

8


⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−

+

+−

=

12

1

2211

2211

22

.

cos..

.cos..

L

gsenL

sengL

v

θ

θθθ

θθθ

&

&&

&&

&

⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−

−

=

02

.

2.

121

11

11 L

L

vC θ

θ

&

&&

&

⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−−−

++−

+−+−

=

12

122

2122

1

221122

22222

1

22112222

12

2

22

.cos2

.2

.

cos..2

.2

.

.cos.cos22

.

LLsenL

gsenLL

senL

sengLsenLL

vC

θθθθθθ

θθθθθθ

θθθθθθθ

&&&&&

&&&&

&&&&&

&

próximo passo é a determinação das forças resultantes no centro de massa de cada

fetuando os cálculos, chegamos a :

Oligamento, segundo a equação abaixo:

11

111 . +

+++

+ = Cii

iii vmF &

E

⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−

−

=

02

212

11

111

11 L

m

Lm

F θ

θ

&

&&

9


⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−−−

++−

+−+−

=

12

1222

21222

12

22211222

222222

12

2221122222

122

22

22

.cos2

.2

.

cos..2

.2

.

.cos.cos22

.

LmLmsenLm

gmsenLmL

msenL

m

sengmLmsenLmLm

F

θθθθθθ

θθθθθθ

θθθθθθθ

&&&&&

&&&&

&&&&&

A seguir, calculamos a somatória dos momentos atuantes nos ligamentos:

).(. 11

11

11

11

11

11

++

++

++

++

++

++ ∧+= i

ii

Cii

ii

ii

Cii

i wIwwIN & Obtemos:

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

=

11

11

.00

zzIN

θ&&

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−

=

22

2221

2221

22

.

..

cos..

zz

zz

zz

I

senI

I

N

θ

θθθ

θθθ

&&

&&

&&

Na seqüência, rearranjamos as equações acima para que fiquem na forma de esforços nas juntas recorrentes do ligamento i+1 para o i :

iii

ii

iii FfRf += +

++111 .

1

1111

11 .. +

++++

++ ∧+∧++= i

iiii

ii

iCi

ii

iiii

ii

i fROFOnRNn Substituindo valores, obtemos:

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−−−

++−

+−+−

==TSR

LmLmsenLm

gmsenLmLmsenLm

sengmLmsenLmLm

Ff

12

1222

21222

12

22211222

222222

12

2221122222

122

22

22

22

.cos2

.2

.

cos..2

.2

.

.cos.cos22

.

θθθθθθ

θθθθθθ

θθθθθθθ

&&&&&

&&&&

&&&&&

10

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional Introduzimos 3 variáveis (R,S,T) correspondendo às 3 linhas da matriz , apenas para simplificar a notação e facilitar a visualização em cálculos futuros.

22 f

⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

++−−

−−−

−++−−−

=

222

222222

1222222

1222

22

12

1222

21222

12

111

22222

321211222

22122

222

11

cos2

)cos(2

)().2(42

.

23cos

2)(

2

222)cos().2(

4cos

2.

θθθθθθθθθθ

θθθθθθ

θθθθθθθθθθ

LmsenLmgmsensenLmsenLm

LmL

msenL

m

LmLmsenLmLmsenLmLm

f

&&&&&

&&&&&

&&&&&&&&&

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

−+−+

−−−−

=FED

senLgmLLmsenLmLmI

senI

LLmLmsenLmI

n

zz

zz

zz

22

2221

122

222

12

22

2222

2221

212122

22

2122

22

122221

22

2.cos

2)2(

84

2cos

44cos

θθθθθθθ

θθθ

θθθθθθθθθ

&&&&&&&

&&

&&&&&&&

Mais 3 variáveis foram introduzidas (D,E,F), também para melhor visualização nos cálculos finais.

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

+−+−−

−−−=

2121

21

12211

212122

11

.cos.4

cos

cos..cos

θθθθθ

θθθθ

senLSLRLmEDsenI

FLSsenLREsenD

n

zz&&

Por fim, o torque na junta i é dado por :

iiT

ii

i zn ).)(=τ De onde concluímos :

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

−+−+

+−+−−=⎥

⎦

⎤⎢⎣

⎡

22

2221

122

222

12

22

2222

2121

21

12211

2

1

2.cos

2)2(

84

.cos.4

cos

θθθθθθθ

θθθθθ

ττ

senL

gmLL

msenL

mL

mI

senLSLRL

mEDsenI

zz

zz

&&&&&&&

&&

11

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional O torque na junta 1 é o resultado da ação externa do motor e dos atritos de Coulomb e viscoso. Na junta 2 há apenas os atritos. Podemos escrever o sistema acima da seguinte forma :

),()(),()( qqFqGqqVqqM &&&& +++=Γ onde : q : vetor de coordenadas generalizadas

)(qM : matriz de massa do manipulador ),( qqV & : matriz de velocidades das juntas

)(qG : matriz de torques gravitacionais ),( qqF & : matriz de torques por atrito

Portanto, o sistema encontrado, escrito na estrutura anterior, é :

444 8444 76

&&

&&

444 3444 21

444444444 8444444444 76

&

&

&

&&

&&

&&

44444444444 344444444444 21

FG

VM

zz

zz

BC

BC

senLgmsenLm

senLLmsenLm

LmILLm

LLmLmLmsenLmI

⎥⎥⎦

⎤

⎢⎢⎣

⎡

+

++

+⎥⎥

⎦

⎤

⎢⎢

⎣

⎡

−+

⎥⎥⎦

⎤

⎢⎢⎣

⎡

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

−

−+

+⎥⎥⎦

⎤

⎢⎢⎣

⎡

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡+++

=⎥⎦

⎤⎢⎣

⎡

2222

1111

22

22

1

12

22

2

2221

222

22

2

2

12

2222

212

221

22

12

21

122

22

211

)sgn(.

)sgn(.

.2

..

0

0)2(8

2)2(

8

4cos

2

cos24

)(4

0

θθ

θθ

θθ

θ

θθ

θθθθ

θ

θ

θ

θθτ onde C1 e C2 são coeficientes de atrito estático, B1 e B2 coeficientes de atrito dinâmico e τ1

o torque do motor. Obs: por ser um sistema sub-atuado, o torque 2τ é igual a zero .

12

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional 2.2 Linearização em torno do equilíbrio instável ( πθ =2 ) 2.2.1 Estratégia 1 : Linearização entrada-saída Por ser o sistema extremamente não linear, ficou no ar a dúvida se seria ou não conveniente uma abordagem de linearização diferente das comumente utilizadas (por Taylor, por exemplo). Pensando nisto, inicialmente foi testada a linearização entrada-saída, que cancela completamente as não-linearidades do sistema. O procedimento é o seguinte: Queremos achar uma relação direta e simples entre a saida e a entrada y 1τ do sistema.

),()(),()( qqFqGqqVqqM &&&& +++=Τ

[ ]FGVTMq −−−= −1&&

[ ]Tq 21 θθ=

[ ] [

[ ] [⎪⎪⎪

⎩

⎪⎪⎪

⎨

⎧

=−−−==

==−−−==

==

−

−

3

124

423

112

211

10

01

xyFGVTMx

xxFGVTMx

xx

θθθθ

&&&

&&

&&&

&&

]

]

Para gerar uma relação direta entre saida e a entrada y 1τ , basta derivar a saida . y

43 xxy == && Como não se relaciona diretamente com a entrada y& 1τ , derivamos novamente e obtemos:

[ ] [ ]FGVTMxy −−−== −14 10&&&

[ ] [ ] [ ]FGVMTMxy ++−== −− 11

4 1010&&&

[ ] [ ] [ FGVMτMxy ++−⎥⎦

⎤⎢⎣

⎡== −− 1

11

4 1001

10&&& ]

13

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional Após algumas expansões, chegamos a:

[ ] [ ] [4444 34444 21

44 344 21

&&&

Nb

FGVMτMxy ++−⎥⎦

⎤⎢⎣

⎡== −− 1

11

4 1001

10 ]

Escolhendo a entrada de controle para ser da forma:

( )Nb −= − ντ 11

onde ν é a nova entrada a ser determinada, chegamos a:

Nby −⋅= 1τ&&

( )[ ] νν =−−⋅= − NNbby 1&& Como se pode observar, a não-linearidade do sistema é cancelada, e obtemos o comportamento de um integrador duplo entre a saida e a nova entrada ν ,

ν=y&& Montando um diagrama de blocos no simulink, temos:

Figura 1 – Esquemático do pêndulo, controlador e linearizador

Uma vez linearizado o sistema, calculamos um controlador que obedecesse às seguintes condições:

• Tempo de subida = 1 s • Sobressinal = 1 %

Como, linearizado, o modelo do pêndulo se comporta como um pólo duplo na

origem (com ganho unitário), a obtenção do controlador ( atraso-avanço) é feita sem maiores dificuldades:

14


452,7309,20)(

+=

sssG

Simulamos o sistema com condição inicial :

[ ] [ ] [ ]TTTxxxx 01.30022114321 == θθθθ &&

e com a referência mantida em 3.15 (próximo a pi). OBS: os parâmetros de atrito empregados na simulação serão determinados em sessões seguintes. Os resultados foram:

0 1 2 3 4 5 6 7 8 9 103.1

3.11

3.12

3.13

3.14

3.15

3.16

3.17Resposta do sistema+controlador linear - x(3)

tempo(s)

desl

ocam

ento

ang

ular

(rad

)

Figura 2 – Resposta θ2 para condições iniciais [ ] [ ]01.3002211 =T

θθθθ &&

15

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional A resposta parece satisfatória, entretanto apenas se não levarmos em conta x(1).

0 1 2 3 4 5 6 7 8 9 10-7

-6

-5

-4

-3

-2

-1

0

1de

sloc

amen

to a

ngul

ar (r

ad)

tempo(s)

Resposta do sistema+controlador linear - x(1)

Figura 3 – Resposta θ1 para condições iniciais [ ] [ ]01.3002211 =T

θθθθ &&

Para manter as condições impostas, o controlador não leva em consideração a posição x(1). Isso acarreta total falta de controle sobre esta variável (que, na realidade, não era o objetivo da linearização). Isto não é interessante, pois o braço horizontal pode girar indefinidamente.

Percebido este fato e sabendo que esta técnica de linearização não é considerada

robusta (ou seja, é muito sensível a erros de modelagem ou imprecisão em parâmetros), decidimos adotar outra abordagem para linearizar o sistema.

16

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional 2.2.2 Estratégia 2 : Linearização por Taylor

Para adoção desta técnica, é conveniente uma mudança de variável ( θ3 ), de forma que o ângulo referente ao deslocamento da barra vertical seja zero quando ela está para cima.

Assim:

23

23

23

θθ

θθ

θπθ

&&&&

&&

−=

−=

−=

Com esta mudança de variável, o modelo fica:

4444 84444 76

&&

&&

4444 34444 21

444444444444 8444444444444 76

&

&

&

&&

&&

&&

4444444444444 34444444444444 21

F

G

VM

zz

zz

bC

bCsenLgm

senLm

senLLmsenLm

LmILLm

LLmLmLmsenLmI

⎥⎥⎦

⎤

⎢⎢⎣

⎡

−+−

++

⎥⎥

⎦

⎤

⎢⎢

⎣

⎡

−−+

+⎥⎥⎦

⎤

⎢⎢⎣

⎡

−⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

−−

−−−−−+

+⎥⎥⎦

⎤

⎢⎢⎣

⎡

−⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

+−

−++−+=⎥

⎦

⎤⎢⎣

⎡

)()sgn(.

)sgn(.)(.

2..

0

0)22(8

))((2

))(22(8

4)cos(

2

)cos(24

)(4

3232

1111

32

2

3

1

13

22

2

3321

233

22

2

3

12

2223

212

321

22

12

21

132

22

21

2

1

θθ

θθθπ

θ

θ

θθπ

θθπθθπ

θ

θ

θπ

θπθπ

ττ

O controle linear deve agir quando a barra vertical está na região de equilíbrio instável (para cima), portanto o sistema acima deve ser linearizado em torno do ponto x :

[ ] [ ] [ ]TTTxxxxx 000033114321 === θθθθ &&

A linearização de um vetor em torno de um ponto é feita por meio de uma aproximação de Taylor :

17


)).(()()( 111 xxxxx −∇+≈ τττr

)).(()()( 222 xxxxx −∇+≈ τττ

r

Efetuados os cálculos, o sistema é simplificado para :

( )

( ) 32322

32

22

21212

2

113212

112

21

21

11

242

24

θθθθτ

θθθτ

&&&&&

&&&&&

bgLmLm

ILLm

bLLmLmLmI

zz

zz

−⎟⎠⎞

⎜⎝⎛−⎟⎟

⎠

⎞⎜⎜⎝

⎛+−⎟

⎠⎞

⎜⎝⎛−=

+⎟⎠⎞

⎜⎝⎛+⎟⎟

⎠

⎞⎜⎜⎝

⎛++=

Para facilitar a visualização, fazemos as seguintes mudanças de variáveis:

2

4

2

4

22

22

22

212

12

21

21

1

gLmD

LmIC

LLmB

LmLmIA

zz

zz

=

+=

=

++=

O sistema anterior com as variáveis substituídas é:

323312

11311

θθθθτ

θθθτ&&&&&

&&&&&

bDCB

bBA

−−−−=

++=

Lembrando que o torque 2τ é igual a zero e passando as equações acima para o domínio de Laplace, obtemos as seguintes funções de transferência:

DsbCsBs

ss

++−

=2

2

2

3

1

)()(

θθ

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

1212

12324

22

1

11 DbsbbADsCbAbsBACs

DsBCssssG

+++++−++

==τθ

18


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

)(12112

223

2

1

32 DbbbADsCbAbsBACs

Bsss

sG+++++−

−==

τθ

Os diagramas de Bode e LGR’s correspondentes às duas últimas funções de transferência são:

LGR de G1(s)

Real Axis

Imag

Axi

s

-2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2

-3

-2

-1

0

1

2

3

Figura 4 – LGR de G1(s)

19


Bode Diagram

Frequency (rad/sec)

Phas

e (d

eg)

Mag

nitu

de (d

B)

-20

0

20

40

60

10-1

100

101

-180

-135

-90

Figura 5 – Diagrama de Bode de G1(s)

LGR de G2(s)

Real Axis

Imag

Axi

s

-8 -6 -4 -2 0 2 4 6 8

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8


20


Bode Diagram

Frequency (rad/sec)

Phas

e (d

eg)

Mag

nitu

de (d

B)

-10

-5

0

5

10

15

10-1

100

101

0

45

90


21


2.3 Linearização em torno do equilíbrio estável ( 02 =θ ) O procedimento para linearização em torno de θ2 = 0 é muito semelhante ao executado para o ponto θ2 = π . Neste caso, entretanto, não faremos a mudança de variável para θ3 , pois o ponto de interesse já é o zero. Utilizando a mesma notação de “A”, “B”, “C” e “D” para simplificação de cálculos, chegamos a :

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

1212

12324

22

1

11 DbsbbADsCbAbsBACs

DsBCssssG

−+−+++−−+

==τθ

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

12112223

2

1

22 DbbbADsCbAbsBACs

BssssG

−+−+++−−

==τθ

Os LGR´s e os Diagramas de Bode correspondentes a estas funções de transferência são:

LGR de G1(s)

Real Axis

Imag

Axi

s

-1.2 -1 -0.8 -0.6 -0.4 -0.2 0

-6

-4

-2

0

2

4

6


22


Bode Diagram

Frequency (rad/sec)

Phas

e (d

eg)

Mag

nitu

de (d

B)

-20

0

20

40

60

10-1

100

101

-180

-135

-90

-45


LGR de G2(s)

Real Axis

Imag

Axi

s

-6 -5 -4 -3 -2 -1 0 1 2 3

-2

-1.5

-1

-0.5

0

0.5

1

1.5

2


23


Bode Diagram

Frequency (rad/sec)

Phas

e (d

eg)

Mag

nitu

de (d

B)

-10

0

10

20

30

40

10-1

100

101

-360

-270

-180

-90


24


3 Componentes Utilizados 3.1 Elétricos O projeto do pêndulo foi inteiramente realizado a partir dos seguintes componentes elétricos :

- 1 motor DC 12V de potência desconhecida; - 1 placa de amplificação de ganho de 2,4V, com ponte H incluída, para

acionamento do motor; - 2 potenciômetros de 3 voltas, de 5kΩ cada, utilizados como sensores de posição

dos ligamentos horizontal e vertical. Posteriormente, verificamos que a precisão de medida angular era de cerca de 2 graus.

No laboratório dispúnhamos de:

- 1 computador equipado com placa de aquisição Lynx, Matlab versão 4.0 e compilador C/C++ Watcom;

- fontes DC de até 25V para alimentação dos potenciômetros e da placa de amplificação.

3.2 Mecânicos O desenho e a construção dos componentes mecânicos foi a etapa mais difícil do projeto, devido ao pouco contato que tivemos com criação de peças mecânicas. A premissa adotada foi a de que as diversas partes do pêndulo pudessem ser desmontáveis, de forma que não houvesse dificuldade caso fosse necessário trocar alguma peça. No Anexo II há o projeto detalhado de cada um dos componentes empregados na montagem. Excetuando-se os rolamentos, que eram de aço, todas as peças construídas (e as polias, utilizadas na redução de 3:1) são de alumínio. A figura a seguir ilustra fielmente o pêndulo invertido rotacional:

25


Figura 12 – Pêndulo invertido rotacional montado O cilindro menor (mais à direita na figura) contém o acoplamento entre o motor e a polia menor. O cilindro maior contém o acoplamento da polia maior e o eixo em que está conectado um dos potenciômetros (o que mede o deslocamento angular da barra horizontal). Na ponta da barra retangular (horizontal), vemos um outro acoplamento. Este contém o potenciômetro que mede o deslocamento angular da barra vertical. Ambos os potenciômetros são idênticos e, como já foi mencionado, têm uma excursão de 1080º. A barra horizontal foi escolhida como sendo um perfil retangular com o intuito de ter seu momento de inércia diminuído. Veremos, na etapa de identificação de parâmetros, que este objetivo não pôde ser alcançado, pois consideramos como sendo o momento de inércia do ligamento horizontal o da própria barra e do acoplamento em sua ponta (cujo momento de inércia é alto).

26

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional No caso da barra vertical, escolhemos um cilindro maciço, para que este tivesse um momento de inércia maior. O comprimento da barra foi estimado em cerca de 30 cm antes da construção. Todos os componentes do pêndulo estão conectados com parafusos pressionando um eixo ao outro, exceto a barra horizontal e seu eixo, que foram prendidos com adesivo industrial. Na foto a seguir vemos o pêndulo invertido rotacional construído e montado:

Figura 13 – Pêndulo invertido rotacional construído e montado

27

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional Experimentamos muitas dificuldades na construção dos componentes. Ainda que com a ajuda do Professor Gilberto e dos técnicos da oficina, algumas peças tiveram que ser refeitas e alteradas diversas vezes. Problemas como eixos que não entravam em furos, rolamentos que travavam, furos fora de centro etc. nos acompanharam por mais de 3 meses. Alguns dos contratempos que tivemos ocorreram devido à nossa pouca habilidade em projeto mecânico, por isso aconselhamos os alunos de engenharia elétrica que se aventurarem em projetar componentes mecânicos:

- a procurarem, como fizemos, a ajuda de um professor de engenharia mecânica; - a pesquisarem como realizar desenhos de fabricação (que é a seqüência de

desenho técnico); - a pesquisarem conexões de eixos (muitas vezes parafusos não são suficientes,

principalmente quando há variações grandes de torque, ou do sentido do torque); - eixos devem ser feitos de aço. Cuidado com o projeto caso optem por eixos de

alumínio. - os encaixes devem ser perfeitos. Um furo com 1 décimo de milímetro fora da

especificação pode se tornar um grande problema. Apesar das dificuldades, conseguimos uma montagem satisfatória. Posteriormente, na fase de implementação, entretanto, verificamos que as imperfeições de construção prejudicaram muito o desempenho, principalmente quanto ao atrito estático (que variava segundo o sentido de rotação e a posição do pêndulo).

28


4 Identificação de Parâmetros 4.1 Parâmetros do Motor A primeira preocupação foi determinar o motor a ser utilizado. Em um primeiro momento, cogitou-se aproveitar o já existente no pêndulo antigo (carrinho), presente em laboratório. Contudo, a retirada deste motor seria desnecessária, pois dispúnhamos de outras alternativas. Foi escolhido um motor DC 12V, fornecido pelo nosso orientador (o professor José Jaime da Cruz). Com o motor em mãos, o primeiro passo foi determinar a potência que ele pode fornecer ao sistema. Isto foi conseguido com a execução de um ensaio no qual a corrente que passa pelo motor e velocidade de rotação do eixo eram determinados a partir de valores de tensão aplicados aos terminais do motor. A resistência da armadura (Ra) foi determinada como sendo 10Ω. Abaixo está ilustrado o modelo do motor adotado:

Raia

e Va

ω

A

Figura 14 – Modelo do Motor

A medição da corrente foi realizada utilizando um multímetro como amperímetro em série com o motor. Outro multímetro serviu como voltímetro. Para medição da velocidade angular no motor utilizamos um estroboscópio.

29


Figura 15 – Ensaio com o estroboscópio – disco girando

Figura 16 – Ensaio com o estroboscópio – disco iluminado

As medições de rotação e de corrente foram feitas variando (de 1 em 1 V) a tensão de 2V a 12V. Para uma maior precisão dos valores medidos, o ensaio foi executado duas vezes: aumentando a tensão a partir de 2V e diminuindo a partir de 12 V.

30


Os resultados obtidos estão ilustrados nas tabelas abaixo: VALORES OBTIDOS VARIANDO-SE A TENSÃO DE 2 A 12 V Tensão

de Entrada

V (V)

Corrente de Entrada

I (mA)

Rotação do motor

(rpm)

Rotação do motor

(rad/s)

Constante de proporcionalidade

da fcem kv(V*s/rad)

Potência do

motor (W)

2 71 175 1099.557429 0.001173199 0.09159 3 72 345 2167.698931 0.001051807 0.16416 4 73.1 515 3235.840433 0.001010248 0.238964 5 72.2 680 4272.566009 0.001001272 0.308872 6 73.2 880 5529.20307 0.000952759 0.385618 7 73.5 1040 6534.512719 0.000958755 0.460478 8 74.3 1220 7665.486075 0.000946711 0.539195 9 74.6 1380 8670.795724 0.000951931 0.615748 10 75.3 1540 9676.105373 0.000955653 0.696299 11 75.4 1720 10807.07873 0.000948082 0.772548 12 75.7 1880 11812.38838 0.000951797 0.851095

média de kv (a partir de 6

V) : 0.000952241

Tabela 2 – Ensaio da característica do motor (aumentando a tensão)

VALORES OBTIDOS VARIANDO-SE A TENSÃO DE 12 A 2 V Tensão

de Entrada V

(V)

Corrente de Entrada I

(mA)

Rotação do motor

(rpm)

Rotação do motor

(rad/s)

Constante de proporcionalidade

da fcem kv(V*s/rad)

Potência do

motor (W)

2 70 175 1099.557429 0.001182294 0.091 3 71 360 2261.946711 0.001012402 0.16259 4 72 520 3267.25636 0.0010039 0.23616 5 72 690 4335.397862 0.000987222 0.30816 6 72 900 5654.866776 0.000933709 0.38016 7 72.5 1060 6660.176426 0.000942167 0.454938 8 74 1220 7665.486075 0.000947102 0.53724 9 74.2 1380 8670.795724 0.000952392 0.612744 10 75.3 1560 9801.769079 0.000943401 0.696299 11 75.4 1720 10807.07873 0.000948082 0.772548 12 76.1 1880 11812.38838 0.000951459 0.855288

média de kv (a partir de 6

V) : 0.000945473

Tabela 3 – Ensaio da característica do motor (diminuindo a tensão)

31

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional Vemos que os dois ensaios forneceram resultados muito parecidos. A partir de 6V de tensão de entrada, influências de ruídos e transitórios passam a ser menos sentidos e os valores relativos a essas tensões mais confiáveis. No gráfico a seguir observamos este fato na curva que representa a constante kv (constante de proporcionalidade da força eletromotriz), que, a partir de 6V, torna-se praticamente constante:

kv(V*s/rad) pela tensão de entrada

0

0.0002

0.0004

0.0006

0.0008

0.001

0.0012

0.0014

2 3 4 5 6 7 8 9 10 11 12

tensão de entrada (V)

kv(V

*s/r

ad)

Serão utilizadas tensões acima de 6V para o acionamento do pêndulo

Figura 17 – Gráfico de kv (constante de proporcionalidade de fcem) pela tensão de entrada

A potência (útil) do motor (presente nas tabelas acima) foi calculada a partir da potência total (do gerador) subtraindo a potência dissipada na resistência de armadura.

32

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional Sua curva em função da tensão de entrada é :

Potência(W) pela tensão de entrada(V)

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

0 2 4 6 8 10 12

tensão de entrada(V)

Pot

ênci

a út

il do

mot

or (W

)

14

Figura 18 – Potência útil do motor em função da tensão em sua entrada

A 8V, a potência útil transmitida ao eixo é de 0,54W, o que é suficiente para acelerar o sistema (cuja massa é aproximadamente 400g) a mais de 1 metro por segundo ao quadrado. Completado o estudo da identificação do motor, passamos à determinação dos demais parâmetros do sistema.

33

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional 4.2 Parâmetros do Sistema 4.2.1 Momentos de Inércia Os primeiros parâmetros a serem determinados foram os momentos de inércia. Calculamos os momentos para o conjunto da barra horizontal (composto pela própria barra e o acoplamento em sua ponta, que chamaremos de acoplamento-barra-vertical) e da barra vertical.

Utilizamos para isso os seguintes dados:

• Massa específica do alumínio: 2700 kg/m3 • Massa da barra horizontal: 64,20 g • Comprimento da barra horizontal: 299,7 cm • Aresta externa da barra: 12,90 mm • Espessura da barra: 1,8 mm

Momento de inércia de uma barra de massa “m”, comprimento “L”, seção quadrada de lado “h” e com o centro de massa transladado de “d” em relação ao eixo de rotação:

222 )(

121 mdhLmI ++=

Tendo conhecimento da espessura da barra, calculamos seu momento de inércia subtraindo do valor correspondente a uma barra maciça (com a aresta e o comprimento indicados acima) o momento de inércia de uma barra menor (de mesmo comprimento L e aresta = 9,30 mm). Concluímos que:

23 .10.007,2 mkgI h

−=

Entretanto, para o ensaio de atrito, montamos o acoplamento-barra-vertical na ponta da barra horizontal para que sua massa ficasse maior e um maior número de oscilações fosse observado. Os dados do acoplamento (que foi aproximado por um cilindro, para efeito de cálculo de momento de inércia) são os seguintes:

• Massa do acoplamento: 170,05 g • Comprimento: 60 mm • Diâmetro: 40 mm

O momento de inércia de um cilindro de massa “m” e comprimento “L” cujo eixo transversal é perpendicular ao eixo de rotação e que tem seu centro de massa deslocado de uma distância “d” do eixo de rotação é :

34


22

121 mdmLI +=

Portanto, o momento de inércia do acoplamento-barra-vertical é

22 .10.238,1 mkgIbv−=

Somando este valor ao momento de inércia da barra horizontal, obtemos o momento de inércia correspondente ao conjunto barra horizontal e acoplamento:

22 .10.439,1 mkgI orizontalconjBarraH−=

Para a barra vertical (cilíndrica), temos os seguintes dados :

• Massa da barra vertical: 102,81 g • Comprimento da barra vertical: 298,2 mm • Diâmetro da barra vertical: 12,1 mm

Utilizando a aproximação de uma barra delgada de comprimento “L” e massa “m”, podemos usar a fórmula abaixo para cálculo do momento de inércia:

2

31 mLI =

Substituindo valores, obtemos:

23 .10.055,3 mkgI calBarraVerti−=

Os demais componentes (eixos, rolamentos, etc) que poderiam contribuir com momentos de inércia foram desconsiderados (isoladamente) por não pertencerem ao modelo. Com este ensaio, pretendemos considerar o sistema como um todo, já levando em conta, portanto, contribuições de rolamentos, eixos etc.

35

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional 4.2.2 Coeficientes de Atrito Viscoso Os ensaios efetuados para a obtenção dos coeficientes de atrito dinâmico consistiram em simular um pêndulo simples no qual o ligamento pendente era o conjunto barra horizontal e acoplamento-barra-vertical (no primeiro ensaio) e somente a barra vertical (no segundo ensaio).

Para melhor precisão dos dados coletados na placa Lynx, aplicamos o seguinte filtro passa-baixa à saída dos potenciômetros:

100100)(+

=s

sG

Este filtro apresenta ganho unitário e atraso de menos de 3 centésimos em relação à resposta sem filtro quan7do um degrau é aplicado em sua entrada. No laboratório, entretanto, verificamos atraso de aproximadamente 1 décimo de segundo.

Os potenciômetros (de 3 voltas) foram alimentados com 10 V, e o gráfico correspondente à resposta livre da barra vertical está mostrado abaixo:

Atrito viscoso - Vertical

-2.500

-2.000

-1.500

-1.000

-0.500

0.0000.00 5.00 10.00 15.00 20.00 25.00

tempo (s)

Tens

ão (V

)

yfiltroysemfiltro

Figura 19 – Comportamento livre da barra vertical – alimentação dos potenciômetros a 10 V

Como se pode perceber, as respostas sem filtro e com filtro são muito semelhantes, de forma que o ruído no sensor não se torna grande problema.

36

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional A envoltória superior da resposta oscilatória anterior é da forma :

TteAtf /.)( −=

onde T é a constante de decaimento procurada.

Para determinação do valor de T, escolhemos 5 períodos da resposta acima e a eles aplicamos o Logaritmo Neperiano , conforme o gráfico a seguir :

Máximos

-1.200

-1.000

-0.800

-0.600

-0.400

-0.200

0.0000 1 2 3 4 5 6 7 8 9 10

tempo (s)

LN d

a Te

nsão

maximoLinear (maximo)

Figura 20 – LN da resposta livre da barra vertical

A derivada da reta presente no gráfico é igual a (-1/T), de onde concluímos :

s83,10T = O período de oscilação para o intervalo considerado é :

sPer 26,1= De posse destes dados, é possível calcular as freqüências naturais amortecida (wn) e não amortecida (wd) :

sradw

wPer d

d

9866,42=→=

π

37


0185,0

9875,4

1

.1

2

=

=∴

⎪⎩

⎪⎨

⎧

−=

−=−

ξ

ξ

ξ

sradw

ww

wT

n

nd

n

Sabendo que a dinâmica do pêndulo pode ser representada por um sistema de segunda ordem do tipo :

ukbI =++ θθθ &&&

que, por sua vez, corresponde à seguinte função de transferência:

22

2

2.)(

nn

n

wswsw

KsG++

=ξ

concluímos:

nwIb ξ2=

radNmsbvertical

410.640,5 −=

Procedimento semelhante foi empregado para determinação do coeficiente de atrito dinâmico da barra horizontal.

As curvas de comportamento livre e os valores calculados a partir delas estão apresentados a seguir:

38


Atrito viscoso - Horizontal

0.000

0.200

0.400

0.600

0.800

1.000

1.200

1.400

1.600

1.800

0.00 2.00 4.00 6.00 8.00 10.00 12.00

tempo (s)

Tens

ão (V

)

yfiltroysemfiltro

Figura 21 - Comportamento livre da barra horizontal – alimentação dos potenciômetros a 10 V

Maximos e Minimos

-2.500

-2.000

-1.500

-1.000

-0.500

0.0000 0.5 1 1.5 2 2.5 3 3.5 4

tempo (s)

LN d

a Te

nsão

maximoLinear (maximo)

Figura 22 – LN da resposta livre da barra horizontal

39


Parâmetros calculados:

sT 68,1=

sPer 04,1=

sradwd 0415,6=

sradwn 0707,6=

0981,0=ξ

radNmsbhorizontal

210.713,1 −=

Verificamos que o amortecimento é muito maior no caso da barra horizontal (o que de fato ocorreu na prática) e, por conseguinte, seu coeficiente de atrito viscoso é muito maior do que o da barra vertical.

40

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional 4.2.3 Coeficientes de Atrito Estático De muita importância para o correto acionamento do motor é o conhecimento do atrito estático presente no sistema composto pelo próprio motor e o eixo conectado à barra horizontal. Como este parâmetro não é modelado nas linearizações, o controle pode ser prejudicado se seu valor for muito alto. Já foi visto que os coeficientes de atrito dinâmico (tanto para a junta da barra horizontal quando a da barra vertical) apresentaram valores relativamente altos, e o que se espera é que o mesmo ocorra com os coeficientes de atrito estático. Foram realizados 2 ensaios, sendo que a barra horizontal era posicionada segundo uma origem arbitrária. Aplicamos na entrada do motor uma rampa de 0V a -3V e outra de 0V a 3V, e o deslocamento angular da barra horizontal foi monitorado. OBS: lembramos que a saída da placa Lynx está conectada a uma placa amplificadora de ganho 2,4V/V. Nos gráficos a seguir, portanto, os valores de tensão devem ser multiplicados por 2,4 para se saber a real voltagem à qual o motor está submetido. Os resultados estão a seguir:

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5-3

-2.5

-2

-1.5

-1

-0.5

0

0.5

1

tempo(s)

Tens

ao n

o P

ot. e

Ten

sao

no M

otor

(V)

Tensao no Pot. da Barra Horizontal e Tensao Aplicada ao Motor

Tensãono Motor

Tensão no Potenc.

Figura 23 – Tensão aplicada no motor (0 a -3V) e correspondente tensão de saída no potenciômetro da barra horizontal.

.

41


0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 50

0.5

1

1.5

2

2.5

3Tensao no Pot. da Barra Horizontal e Tensao Aplicada ao Motor

tempo(s)

Tens

ao n

o P

ot. e

Ten

sao

no M

otor

(V)

Tensãono Motor

Tensão no Potenc.

Figura 24 – Tensão aplicada no motor (0 a 3V) e correspondente tensão de saída no potenciômetro da barra horizontal.

As curvas anteriores revelam uma característica indesejada do sistema montado: o atrito estático não é constante em toda a excursão do braço horizontal. Quando aplicada tensão negativa, o atrito observado é muito maior do que quando a tensão é positiva. Aplicada tensão positiva, há valores diferentes de atrito conforme a posição angular. Conforme citado anteriormente, a construção não exata das peças (ainda que com imperfeições muito pequenas) resultou em uma montagem que apresenta problemas, como este. Para sanar o problema da diferença de atrito para rotações de sentidos diferentes, implementamos o controle de forma que uma tensão correspondente ao atrito de destaque fosse aplicada conforme o sentido de rotação. Este procedimento não surtiu o efeito esperado, pois a tensão de entrada no motor passou a ser muito alta. Decidimos desconsiderar o atrito estático para a implementação do algoritmo. Esta decisão foi baseada no fato de que, durante a execução do controle, o sistema do pêndulo vibra muito, o que (sabe-se) reduz o efeito do atrito estático.

42


5 Controle 5.1 Elevação A estratégia de elevação foi executada levando em consideração o deslocamento e velocidade dos ligamentos, de forma a sempre aumentar a energia do sistema. No Anexo I encontra-se o algoritmo utilizado (motion.m), e seu funcionamento está descrito abaixo (considerando o sistema de coordenadas adotado): - se a posição da barra vertical for positiva e sua velocidade também for positiva, significa que o pêndulo está subindo, portanto ele é freado (torque negativo é aplicado); - se a posição da barra vertical for positiva e sua velocidade for negativa, significa que o pêndulo está descendo, portanto ele é acelerado (torque positivo é aplicado); - se a posição da barra vertical for negativa e sua velocidade também for negativa, significa que o pêndulo está subindo, portanto ele é freado (torque negativo é aplicado); - se a posição da barra vertical for negativa e sua velocidade for positiva, significa que o pêndulo está descendo, portanto ele é acelerado (torque positivo é aplicado); O torque a ser aplicado (negativo ou positivo) varia conforme o deslocamento angular da barra vertical, sendo maior quando este está distante de π (ou de -π) e menor conforme o pêndulo se aproxima da posição de equilíbrio instável. Há ainda uma imposição para a excursão da barra horizontal, de forma que ela não gire indefinidamente. No algoritmo, esta característica está destacada no último trecho, onde a amplitude de deslocamento da barra horizontal é limitada em 30 graus. Por sintonia deste parâmetro, podemos limitar o deslocamento real (que no caso da simulação foi de cerca de 60 graus).

43

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional O esquemático de Simulink utilizado para simulação foi:

Toque inicial

entrada

Enable

x4

x3 = y2

x2

x1 = y1

x1

u

x3f ator

memória x3

u

teta1

teta2

x3

x1

u

x1-out

StandUP

Saturation

estados entrada

entrada u

estados saída

PÊNDULO

em

Figura 25 – Esquemático utilizado na simulação da elevação

No diagrama acima foi empregada uma saturação de 0,5N.m no torque do motor. O bloco “PÊNDULO” representa o modelo não-linear, o “memória x3” guarda a última iteração do estado x3 (velocidade da barra vertical) e aumenta o fator do torque a ser aplicado (a cada vez que o pêndulo inicia uma queda) e, finalmente, o bloco “Stand UP” contém a rotina motion.m, que melhor seleciona o valor do torque. O resultado da simulação está ilustrado a seguir:

44


0 0.5 1 1.5 2 2.5 3-1.5

-1

-0.5

0

0.5

1θ1 e θ2 durante elevacao

θ1

0 0.5 1 1.5 2 2.5 3-2

-1

0

1

2

3

4

tempo (s)

Pos

icao

Ang

ular

(rad

)

θ2

Figura 26 – Respostas dos deslocamentos angulares dos ligamentos durante a elevação Observamos que as respostas confirmam o que era esperado e, em menos de 3 segundos o pêndulo é invertido. O deslocamento do braço horizontal está claramente limitado em cerca de 1 radiano. O modelo de pêndulo utilizado foi o não-linear, portanto espera-se que as respostas sejam as mais fiéis possíveis. A elevação não foi implementada em laboratório, conforme será comentado no item “Conclusões”, mas o sucesso da simulação mostra que o método desenvolvido pode ser empregado em ensaios reais.

45

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional 5.2 Equilíbrio 5.2.1 Controle Linear em Tempo Contínuo Para o controle de equilíbrio, tomamos como aproximação do sistema as equações derivadas no item 2.2.2 deste relatório. Devido à proximidade de pólos e zeros no semi-plano direito presentes nas funções de transferência G1(s) e G2(s), seus cancelamentos através de técnicas clássicas de controle são arriscados, pois se o cancelamento não for perfeito haverá mais pólos e zeros no lado direito (aumentando a instabilidade do sistema). Portanto, inicialmente, adotamos uma estratégia de controle baseada na realimentação de estados. Para um sistema linear invariante no tempo

BuAxx +=& uma lei de controle representando uma realimentação de estados tem a forma :

Kxxkxkxkxku nn −=++++−= )( 332211 K Substituindo na equação anterior:

xAxBKABKxAxx c=−=−= )(& Os pólos do sistema são definidos por:

0=− cAIλ Teoricamente, estes pólos podem ser quaisquer. Uma maneira adequada de determinar uma resposta rápida e estável do sistema é utilizando a técnica numérica LQR(Linear Quadratic Regulator). O Matlab dispõe deste recurso (função “lqr”).

Para executar este recurso, entretanto, é necessário que o sistema composto por G1(s) e G2(s) seja expresso no espaço de estados. Isto é feito identificando os estados e nas duas equações abaixo:

2x& 4x&

323312

11311

θθθθτ

θθθτ&&&&&

&&&&&

bDCB

bBA

−−−−=

++=

Assim, chegamos a:

46


1

2

2

4

3

2

1

22

221

22

221

4

3

2

1

)(

0)(

0

)()()(0

1000)()()(

00010

τ

⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−

−−

+

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

−−−−

−−

−−

−=

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

ACBB

ACBC

xxxx

ACBAb

ACBAD

ACBBb

ACBBb

ACBBD

ACBCb

xxxx

&

&

&

&

Entretanto, na implementação real utilizamos como entrada valores de tensão, e não torque. Para efetuar a mudança de variável, sabemos que:

⎪⎩

⎪⎨⎧

=

+=

IkkkkRIV

gv

gv

1

1.τ

θ&

onde: I = corrente de armadura do motor R = resistência de armadura do motor kv = constante de força eletromotriz do motor kg = redução do motor para o eixo da barra horizontal Sabendo que kg = 3 (razão entre as polias utilizadas), o sistema é modificado para:

V

ACBBk

ACBCk

xxxx

ACBAb

ACBAD

ACBR

BkBb

ACBBb

ACBBD

ACBR

CkCb

xxxx

v

v

v

v

⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−

−−

+

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−−−

+−

−−

−−

−

−

=

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

)(3

0)(

30

)()()(

9

0

1000)()()(

9

0

0010

2

2

4

3

2

1

22

22

1

22

22

2

1

4

3

2

1

&

&

&

&

A representação das saídas é:

⎥⎦

⎤⎢⎣

⎡+

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

⎥⎦

⎤⎢⎣

⎡=⎥

⎦

⎤⎢⎣

⎡00

01000001

4

3

2

1

2

1

xxxx

yy

Conhecendo os parâmetros do sistema, a matriz de ganhos fornecida pela rotina “lqr” do Matlab é :

47


[ ]5647.02615.11149.00010.01000=K Além do controlador, é necessário um observador de estados, pois no sistema real só dispomos dos estados x1 e x3 . O observador é calculado de forma semelhante ao controlador (também pelo LQR). A injeção de saída, portanto, é :

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

−−

−−

=

3109,92159,06468,41701,09346,07555,01701,01171,2

L

O diagrama a seguir ilustra o sistema linearizado, contendo controlador e observador. Uma perturbação de amplitude 0.2 radiano e 0.2 segundo foi adicionada à saída y2 :

48


Figura 27 – Diagrama de blocos do sistema linearizado + controlador + observador de estados

A seguir, as respostas simuladas dos estados e dos estados estimados:

0 2 4 6 8 10 12 14 16 18 20-0.8

-0.6

-0.4

-0.2

0

0.2Estado x1 Real e Estimado

med

ida

de d

eslo

cam

ento

ang

ular

(rad

)

RealEstimado

0 2 4 6 8 10 12 14 16 18 20-5

0

5

10

15

20x 10

-4

tempo (s)

Erro de Estimacao

Figura 28 – Respostas de x1 para o sistema linearizado + controlador + observador de estados

49


0 2 4 6 8 10 12 14 16 18 20-0.06

-0.04

-0.02

0

0.02

0.04m

edid

a de

des

loca

men

to a

ngul

ar (r

ad)

Estado x3 Real e Estimado

tempo (s)0 2 4 6 8 10 12 14 16 18 20

-0.06

-0.04

-0.02

0

0.02

RealEstimado

Erro de Estimacao

Figura 29 – Respostas de x3 para o sistema linearizado + controlador + observador de estados Vemos que o estimador funciona adequadamente para deslocamentos pequenos. Assim, é possível implementar o controle com realimentação de estados. No laboratório, entretanto, sabemos que o Matlab não é executado em tempo real e que, como o sistema necessita de atuações rápidas e precisas, o desempenho pode ser prejudicado. A perda de performance foi verificada nos ensaios (era visível o atraso na resposta). Decidimos, pois, implementar toda a lógica de controle em C (utilizando o compilador Watcom), de modo que a execução se desse em tempo real.

50

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional 5.2.2 Controle Linear em Tempo Discreto O controle em tempo discreto foi realizado utilizando LQR para tempo discreto. A partir da representação do sistema em espaço de estados no tempo contínuo, obteve-se a descrição discreta do sistema em espaço de estados, com segurador de ordem zero e com um tempo de amostragem T de 0,01 segundos, ou seja, uma freqüência de amostragem de 100Hz.

Para obter o modelo discreto utilizamos a rotina c2d.m do Matlab

sistemaC=ss(A,B,C,D) sistemaD=c2d(sistemaC,0.01,'zoh') [phi,gam,H,J]=ssdata(sistemaD)

Calculamos então a matriz de realimentação de estados KD para o caso discreto, através das rotinas de LQR discreto do Matlab: KD=dlqr(phi,gam,eye(4),1)

[ ]558.74-1248.2-113.7-0.9780-=KD A lei de controle discreta fica da forma:

( ))()()()()( 44332211 kxkkxkkxkkxkku +++−= Como uma primeira tentativa, não utilizamos um observador para obter os estados a que não temos acesso. Ao invés disso, eles foram aproximados através de equações de diferenças, pois, e são, respectivamente, e . 2x 4x 1x& 3x&

Tkxkxkxkx )1()()()( 11

12−−

≈= & e T

kxkxkxkx )1()()()( 3334

−−≈= &

Ao testarmos na planta real, percebemos que esta não foi uma boa solução, pois o ruído dos potenciômetros era alto, o que fazia com que as derivadas discretas resultassem em valores cerca de 100 vezes maiores que o esperado. Decidimos então incorporar um observador no nosso sistema de controle discreto. O observador é calculado de forma similar a do controlador. Calculamos então a matriz de injeção de saída LD:

LD= dlqr(phi',H',eye(4),eye(2))'

51


⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

=

0.64680.08650.62810.00140.07400.31610.00130.6232

LD

O observador, estimador, que obtemos é do tipo preditivo, ou seja, estimamos os estados do instante seguinte. A equação do estimador preditivo é:

[ ])()()()()1( kxHkyLDkukxkx −+Γ+Φ=+

A lei de controle, agora com o estimador de estados fica:

( ))()()()()( 44332211 kxkkxkkxkkxkku +++−=

52


6 Implementação e Testes O sistema com o controlador discretizado foi executado diversas vezes. Em todas, a estimação de estados apresentou resultados muito bons (como será mostrado a seguir), com erros praticamente nulos em relação aos estados possíveis de serem medidos. Entretanto, os problemas relacionados às deficiências na montagem vieram à tona. Realmente, apenas uma pequena região de deslocamento angular próxima à origem do sistema de coordenadas da barra horizontal apresentava atrito suficientemente uniforme para o equilíbrio do pêndulo. Fora desta região, o pêndulo dificilmente se equilibrava por mais de 5 segundos, pois o torque do motor era fraco demais para romper o atrito ou então, tentando corrigir este problema, aplicava um torque muito forte, o que impossibilitava a correção no momento seguinte. Os gráficos a seguir ilustram o seguinte experimento : desconectada a atuação do motor do pêndulo, manualmente oscilamos o pêndulo invertido em torno da origem (lembrando: a variável que estamos utilizando para descrever o deslocamento angular da barra vertical é 3θ ) e verificamos o desempenho do estimador. Os resultados estão ilustrados abaixo: OBS: Os valores REAIS de tensão aos quais o motor está submetido são os apresentados nos gráficos multiplicados por 2,4, pois o esforço de controle já foi dividido por este valor no próprio algoritmo.

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5-5

-4

-3

-2

-1

0

1

2

3

4

5Estado X3 Real e Estimado

tempo (s)

Pos

.Ang

.(rad

) e T

ensa

o(V

)

RealEstimadoEsforco de Controle

Figura 30– Estado x3 real, estimado e esforço de controle quando o motor está desconectado

53


0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5-5

-4

-3

-2

-1

0

1

2

3

4

5

tempo (s)

Pos

.Ang

.(rad

) e T

ensa

o(V

)

Estado X1 Real e Estimado


Figura 31 – Estado x1 real, estimado e esforço de controle quando o motor está desconectado

Notamos que os estados estimados são praticamente idênticos aos reais, e que o esforço de controle age no sentido de corrigir o erro em x3. Naturalmente, pelo fato do motor estar desconectado, o valor do esforço aumenta até a saturação da placa Lynx (cerca de 5 V e -5V). Este ensaio nos mostrou outro problema: quando os estados estão com valores nulos, o valor do esforço de controle também deveria estar próximo a zero. O que os gráficos mostram são valores altos de tensão (maiores do que a tensão necessária para vencer o atrito estático, de cerca de 3V) na região em que os estados são zero. Verificamos que há uma imprecisão muito alta com a linearidade dos potenciômetros. A calibração que empregamos orientava a origem para baixo, conforme 2θ (lembrando :

3θ = π- 2θ , portanto medíamos 2θ e no algoritmo subtraíamos de π este valor para realimentarmos 3θ ), de forma que quando a barra estivesse para cima a medida fornecesse 180 graus. Na realidade, o valor que estava sendo medido era cerca de 169 graus, portanto com um erro de 11 graus. Este erro é inaceitável, pois o algoritmo entende que a barra está caindo para um dos lados. Resolvemos este problema orientando a origem para cima, ou seja, já efetuando a medida direta de 3θ .

54

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional Os gráficos a seguir ilustram os resultados obtidos:

0 1 2 3 4 5 6 7 8-5

-4

-3

-2

-1

0

1

2

3

4

5P

os.A

ng.(r

ad) e

Ten

sao(

V)

tempo (s)

Estado X3 Real e Estimado e Esforco de Controle


Figura 32 – Estado x3 real, estimado e esforço de controle

0 1 2 3 4 5 6 7 8-5

-4

-3

-2

-1

0

1

2

3

4

5Estado X1 Real e Estimado e Esforco de Controle

tempo (s)

Pos

.Ang

.(rad

) e T

ensa

o(V

)



55


Observamos que o pêndulo, abandonado do equilíbrio instável (para cima), permanece equilibrado. O esforço de controle é alto mesmo nas condições mostradas nos gráficos, pois para vencer o atrito é necessária uma tensão de mais de 3V(multiplicada por 2,4). O resultado é uma vibração incessante do pêndulo. O ensaio anterior foi conseguido após diversas tentativas, devido à dificuldade do pêndulo em se manter elevado em regiões em que o atrito na junta do braço horizontal fosse grande ou variasse muito. Para tentar melhorar o desempenho do sistema, recebemos orientação para empregar uma barra vertical de comprimento maior. Dispúnhamos de uma segunda barra vertical, de dimensões idênticas à primeira, de forma que acoplamos a ponta de uma em outra e obtivemos uma barra vertical de comprimento e massa duas vezes maiores do que as da original. Recalculando todas as matrizes necessárias para o controle do pêndulo na nova condição, implementamos o algoritmo. As curvas a seguir exibem os resultados do ensaio com a barra maior:

0 1 2 3 4 5 6 7 8-5

-4

-3

-2

-1

0

1

2

3

4

5Estado X3 Real e Estimado e Esforco de Controle

tempo (s)

Pos

.Ang

.(rad

) e T

ensa

o(V

)



56


0 1 2 3 4 5 6 7 8-5

-4

-3

-2

-1

0

1

2

3

4

5

Pos

.Ang

.(rad

) e T

ensa

o(V

)

tempo (s)

Estado X1 Real e Estimado e Esforco de Controle



Com os ensaios percebemos que há uma estabilidade maior, embora ainda haja vibração da barra. O esforço de controle apresentou valores muito altos, que serão objeto de discussão nas conclusões.

57


7 Conclusões O projeto envolveu diversas etapas diferentes, podendo ser considerado inclusive como dois grandes projetos independentes: a construção do aparato e o controle do pêndulo. A fase de construção mecânica ofereceu a maior dificuldade, pois os componentes não eram triviais e de fácil elaboração. Os problemas envolvidos no torneamento das peças também consumiram muito tempo. Devido ao alto custo dos encoders, que serviriam como sensores de posição, optou-se pela utilização de potenciômetros de precisão, embora apresentassem ruído considerável e não fossem plenamente lineares. O motor também ofereceu problemas na implementação, pois muitas vezes não fornecia o torque necessário para correção de movimentos mais bruscos da barra vertical. Suspeitamos que a identificação dos parâmetros possa conter erros. Isto não pôde ser verificado, pois o motor está permanentemente ligado ao acoplamento que o sustenta. Este último problema nos fez refazer o controle diversas vezes, aumentando e diminuindo o ganho da tensão que alimentaria o motor. Em alguns casos, conseguimos um equilíbrio satisfatório da barra vertical. A implementação da elevação não foi possível porque demos prioridade ao desenvolvimento e testes do controle de equilíbrio e também porque (e que reforçou a nossa decisão) suspeitamos que o motor não conseguiria fornecer o torque necessário para o procedimento de elevação. Devemos ressaltar novamente que o melhor desempenho de equilíbrio foi conseguido quando dobramos o tamanho da barra. Finalmente, apesar dos inúmeros contratempos, devemos dizer que este trabalho nos foi extremamente útil, pois enfim nos deparamos com um projeto que englobasse diversas áreas de engenharia e, sobretudo, envolvesse um uso extenso das técnicas de controle aprendidas durante o curso.

58


8 Cronogama Cronograma cumprido.

1. Revisão bibliográfica. 2. Modelagem matemática da dinâmica do pêndulo. 3. Simulação dos modelos obtidos no Matlab/Simulink. 4. Reconstrução mecânica do sistema. 5. Cálculo e identificação dos parâmetros. 6. Projeto do(s) controlador(es). 7. Implementação e testes. 8. Documentação: conclusões e redação final.

Cronograma

Etapas Mar Abr Mai Jun Jul Ago Set Out Nov Dez 1 2 3 4 5 6 7 8

Realizado Não Realizado Continuado

59


9 Bibliografia [1] Ogata, K. “Modern Control Engineering”, 3rd Edition, Prentice-Hall, 1997. [2] Duarte, Cassius R., Zanata, Diogo R. P. “Controle de Pêndulo Invertido”. Escola

Politécnica da USP. [3] Stimac, Andrew K. “Standup and Stabilization of the Inverted Pendulum”

Massachusetts Institute of Technology, 1999. [4] Potsaid, B., Wen, J.T. “Optimal mechanical design of a rotary inverted pendulum.”

International Robotics and Systems (IRoS) Conference, Lausanne, Switzerland, Sept 2002.

[5] K. Yoshida, “Swing-up control of an inverted pendulum by energy-based methods”,

Proceedings of the American Control Conference 1999, pp. 4045-4047. [6] Franklin, G.F., Powell, J.D., Workman, M. “Digital Control of Dynamic Systems”,

3rd Edition, Addison Wesley, 1998.

60


10 Anexo I (Scripts de Matlab e Códigos Fonte) 10.1 Modelo Não Linear do Pêndulo Invertido (sistemaEuler.m) function[xp]=sistema_2(z) %obtencao dos valores dos estados n=4; m=1; for I=1:n, x(I)=z(I); end for I=(n+1):(n+m), u(I-n)=z(I); end z/pi*180; %parametros m1 = 0.23425; m2 = 0.10281; L1 = 0.3; L2 = 0.3; B1 = 0.03007; B2 = 0.000748; C1 = 0.04; C2 = 0.001; Izz1 = 0.01439; Ixx2 = 0; Iyy2 = 0; Izz2 = 0.003055; g = -9.78; kv=0.00095; R=10; %calculo de M M=[Izz1+m2*L2^2*(sin(x(3)))^2/4+m1*L1^2/4+m2*L1^2 m2*L2*L1*cos(x(3))/2; m2*L1*L2*cos(x(3))/2 Izz2+m2*L2^2/4]; T=[(3*kv*u(1))/R-9*(kv^2)*x(2)/R; 0 ]; M_inv=inv(M); V=[ m2*L2^2*sin(2*x(3))*x(4)*x(2)/8-m2*L1*L2*sin(x(3))*x(4)^2/2; -m2*L2^2*(sin(2*x(3))*x(2)^2)/8]; %ja corrigido de L1 para L2 G=[ 0; -m2*g*L2*sin(x(3))/2];

61

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional F=[ C1*sign(x(2))+B1*x(2); C2*sign(x(4))+B2*x(4)]; Result=M_inv*[T-V-G-F]; %estados finais xp(1) = x(2); xp(2) = [1 0]*Result; xp(3) = x(4); xp(4) = [0 1]*Result; 10.2 Modelo Linear: Funções de Transferência e Representação em

Espaço de Estados (snormal_V.m) %parametros m1 = 0.23425; m2 = 0.10281; L1 = 0.3; L2 = 0.3; B1 = 0.01713; B2 = 0.000564; C1 = 0.04; C2 = 0.001; Izz1 = 0.01439; Ixx2 = 0; Iyy2 = 0; Izz2 = 0.003055; g = -9.78; kv=0.00095; Ra=10; %simplificaçoes de termos A=Izz1 + m1*(L1^2)/4 + m2*L1^2; B=m2*L1*L2*0.5; C=Izz2 + m2*(L2^2)/4; D=m2*g*(L2^2)/4; %TF TT1 numTT1 = [ C B2 D]; denTT1 = [ A*C-B^2 A*B2+C*B1 A*D+B1*B2 B1*D 0]; TT1=tf(numTT1,denTT1) %TF TT2 numTT2 = [ -B 0 ]; denTT2 = [ A*C-B^2 A*B2+C*B1 A*D+B1*B2 B1*D]; TT2=tf(numTT2,denTT2) %TF XSI numXSI = [ C-B*(L2/L1)*0.5 B2 D]; denXSI = [ A*C-B^2 A*B2+C*B1 A*D+B1*B2 B1*D 0]; XSI=tf(numXSI,denXSI); %passa o sistema para a representacao em espaco de estados AA=[0 (B^2-A*C) 0 0; 0 C*B1-(9*kv^2*C)/Ra -B*D -B*B2;

62

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional 0 0 0 (B^2-A*C); 0 -B1*B+9*kv*B/Ra A*D A*B2]; AA=AA/(B^2-A*C); BB=[0; -3*kv*C/Ra ; 0 ; 3*kv*B/Ra]; BB=BB/(B^2-A*C); CC=[1 0 0 0; 0 0 1 0]; DD=[0;0]; %calcula a controlabilidade do sistema CO=CTRB(AA,BB); rank(CO) svd(CO) %calcula o controlador pelo LQR Q=eye(4); R=1; [K,S,E] = LQR(AA,BB,Q,R); K=-K %calcula o observador pelo LQR Q=3*eye(4); R=eye(2); [L,S,E] = LQR(AA',CC',Q,R); L=L' %discreto sysC=ss(AA,BB,CC,DD); sysD=c2d(sysC,0.01,'zoh'); [phi,gam,H,J]=ssdata(sysD); KD=dlqr(phi,gam,eye(4),1); LD=dlqr(phi',H',eye(4),eye(2)); LD=LD'; syms x1 x2 x3 x4 y1 y2 u; X=[x1;x2;x3;x4]; Y=[y1;y2]; Xp=phi*X+gam*u+LD*(Y-H*X); 10.3 Elevação do Pêndulo (motion.m) function[excitacao]=motion(z) %obtencao dos valores dos estados

63

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional n=2; for I=1:n, x(I)=z(I); end x(3)=z(3); %mudança de nomenclatura: %dentro desta funçao o estado x(3) e chamado "x(1)", o estado "x(4)" e %chamado "x(2)" e x(1) de "x(3)" if(x(1)>0)%pendulo esta no quadrante positivo if(x(2)>=0) %pendulo esta subindo u=-0.5*abs(pi-abs(x(3))); %breca end if(x(2)<0) %pendulo esta descendo u=1*abs(pi-abs(x(3))); %acelera end end if(x(1)<=0)%pendulo esta no quadrante negativo if(x(2)<=0) %pendulo esta subindo u=0.5*abs(pi-abs(x(3))); %breca end if(x(2)>0) %pendulo esta descendo u=-1*abs(pi-abs(x(3))); %acelera end end %limite do braco horizontal if (x(3)<=-pi/6) u=1.4*abs(pi-abs(x(3))); end if(x(3)>=pi/6) u=-1.4*abs(pi-abs(x(3))); end excitacao(1)=u; 10.4 Código Fonte do Controlador em C (lqrdig.c) #include <stdio.h> // funcoes de entrada e saida #include <conio.h> // usado apenas para a funcao "getch()" #include "labaut.c" //inclui a biblioteca com as rotinas da placa #define NUMERO_DE_AMOSTRAS 1000 // numero de amostras void main(void)

64

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional FILE *arq_c0, *arq_c1, *arq_c2, *arq_c5, *arq_c6, *arq_c7, *arq_c8; int i; // contador. float valor[16]; // valor instantaneo dos dados coletados float X1[NUMERO_DE_AMOSTRAS]; float X3[NUMERO_DE_AMOSTRAS]; float U[NUMERO_DE_AMOSTRAS]; float E1[NUMERO_DE_AMOSTRAS]; float E3[NUMERO_DE_AMOSTRAS]; float E2[NUMERO_DE_AMOSTRAS]; float E4[NUMERO_DE_AMOSTRAS]; float R; float K1,K2,K3,K4; float offset1,offset3; float Kangulo; /* A rotina "setclock" inicializa a placa e programa a frequencia de amostragem * * a ser utilizada. Esta rotina deve ser invocada uma unica vez no inicio de cada * * programa. */ setclock(100); /* reseta a placa */ readA2D(valor); writeD2A(0.0,0.0); readA2D(valor); printf("\n Offset(Volts) Horizontal: %f",(-valor[1])); printf(" Vertical: %f",(-valor[3])); /* calcula offset */ offset1=-valor[1]; offset3=-valor[3]; printf("\n Controlador com Realimentacao de Estados\n"); printf("Pressione qualquer tecla para iniciar a coleta de dados.\n"); getch(); // aguarda que uma tecla seja pressionada printf("Coletando dados...\n"); /* incicializa variaveis */ X1[0]=0;/* Y */ X3[0]=0;/* Y'*/ // X2[0]=0; // X4[0]=0; U[0]=0;

65

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional R=0; E1[0]=0;/* Y */ E3[0]=0;/* Y'*/ E2[0]=0; E4[0]=0; K1=-0.97897; K2=-113.7; K3=-1248.2; K4=-558.74; /* Ganho [rad/volt] */ Kangulo=6.2832/(2.0508-(-0.9473)); for (i=0;i<=NUMERO_DE_AMOSTRAS;i++) readA2D(valor); /* orientacao dos pot*/ X1[i]=Kangulo*(-valor[1]-offset1); X3[i]=-Kangulo*(-valor[3]-offset3); /* esforco de controle */ U[i]= (-K1*E1[i] -K2*E2[i] -K3*E3[i] -K4*E4[i]); if(U[i]>4.9) U[i]=4.9; if(U[i]<-4.9) U[i]=-4.9; /* escreve os dados */ writeD2A(U[i], U[i]); /* Estimador de Estados Preditivo */ E1[i+1]= 0.3768*E1[i] +0.0010*E2[i] -0.0013*E3[i] +0*E4[i] +0*U[i] +0.6232*X1[i] +0.0013*X3[i]; E2[i+1]= 0.9932*E2[i] -0.0818*E3[i] +0.0002*E4[i] +0.0001*U[i] +0.3161*X1[i] -0.3161*E1[i] +0.0740*X3[i]; E3[i+1]= 0*E2[i] +0.3721*E3[i] +0.0100*E4[i] +0*U[i]

66

PTC2531 Elevação e Estabilização Vertical do Pêndulo Invertido Rotacional +0.0014*X1[i] -0.0010*E1[i] +0.6281*X3[i]; E4[i+1]= 0.0056*E2[i] -0.5980*E3[i] +0.9990*E4[i] -0*U[i] +0.0865*X1[i] -0.0865*E1[i] +0.6468*X3[i]; /* reseta a placa */ writeD2A(0.0,0.0); printf("Salvando arquivos...\n"); arq_c0=fopen("lqrx1.dat","w"); arq_c1=fopen("lqru.dat","w"); arq_c2=fopen("lqrx3.dat","w"); arq_c5=fopen("lqre1.dat","w"); arq_c6=fopen("lqre2.dat","w"); arq_c7=fopen("lqre3.dat","w"); arq_c8=fopen("lqre4.dat","w"); for (i=0;i<NUMERO_DE_AMOSTRAS;i++) fprintf(arq_c0,"%g\n",X1[i]); fprintf(arq_c1,"%g\n",U[i]); fprintf(arq_c2,"%g\n",X3[i]); fprintf(arq_c5,"%g\n",E1[i]); fprintf(arq_c6,"%g\n",E2[i]); fprintf(arq_c7,"%g\n",E3[i]); fprintf(arq_c8,"%g\n",E4[i]); fclose(arq_c0); fclose(arq_c1); fclose(arq_c2); fclose(arq_c5); fclose(arq_c6); fclose(arq_c7); fclose(arq_c8); printf("Fim.\n");

67


11 Anexo II (Desenho das Peças)

Desenho das Peças

68

Documents

“Elevação e Estabilização Vertical do Pêndulo Invertido ...dulley/dulley_docs/Pendulo_Invertido_Rotacional... · 2.1 Modelagem Dinâmica ... Além da análise de modelagem