UNIVERSIDADE FEDERAL FLUMINENSE CENTRO TECNOLÓGICO

ESCOLA DE ENGENHARIA PÓS-GRADUAÇÃO EM ENGENHARIA MECÂNICA

LEONARDO RANGEL DIAS FERREIRA

TÍTULO: DESENVOLVIMENTO DE UMA BANCADA DE TESTES PARA VALIDAÇÃO DE UM TANQUE DE ONDAS NUMÉRICO

NITERÓI 2010

LEONARDO RANGEL DIAS FERREIRA

TÍTULO: DESENVOLVIMENTO DE UMA BANCADA DE TESTES PARA VALIDAÇÃO DE UM TANQUE DE ONDAS NUMÉRICO

Dissertação apresentada ao Curso de Pós-Graduação em Engenharia Mecânica da Universidade Federal Fluminense, como requisito parcial para a obtenção do Grau de Mestre. Área de Concentração: Mecânica dos Fluidos.

Orientadores: Prof. JOSÉ ANDRÉS SANTISTEBAN LARREA, D.Sc. Prof. ROGER MATSUMOTO MOREIRA, Ph.D.

NITERÓI 2010

LEONARDO RANGEL DIAS FERREIRA

TÍTULO: DESENVOLVIMENTO DE UMA BANCADA DE TESTES PARA VALIDAÇÃO DE UM TANQUE DE ONDAS NUMÉRICO

Dissertação apresentada ao Curso de Pós-Graduação em Engenharia Mecânica da Universidade Federal Fluminense, como requisito parcial para a obtenção do Grau de Mestre. Área de Concentração: Mecânica dos Fluidos.

BANCA EXAMINADORA

Prof. JOSÉ ANDRÉS SANTISTEBAN LARREA, D.Sc. – Orientador Universidade Federal Fluminense

Prof. ROGER MATSUMOTO MOREIRA, Ph.D. – Orientador Universidade Federal Fluminense

Prof. ANTÔNIO LOPES GAMA, D.Sc. Universidade Federal Fluminense

Prof. GUSTAVO CÉSAR RACHID BODSTEIN, Ph.D.

Universidade Federal do Rio de Janeiro

Prof. JULIO TOMÁS AQUIJE CHACALTANA, D.Sc. Universidade Federal do Espírito Santo

NITERÓI 2010

RESUMO

Em reservatórios que transportam e armazenam os mais diversos tipos de líquidos, ondas com grande energia podem ser criadas a partir de uma simples excitação do sistema. Quando embarcados em um veículo, esses reservatórios poderão eventualmente sofrer danos pelas ondas e até poderão dificultar a manobrabilidade desse veículo. Em particular, o “Sloshing” é um fenômeno onde uma onda com grande quantidade de energia concentrada impacta sobre uma parede de um reservatório e é por isso um dos fenômenos cuja simulação numérica é muito apreciada. No laboratório de dinâmica de fluidos computacional da UFF, este tipo de simulações vem sendo efetuado com sucesso utilizando o método de volumes finitos, sendo que sua validação foi executada, até o momento, comparando com resultados experimentais obtidos em outras partes do mundo. Com isto em mente, o presente trabalho teve como objetivo o desenvolvimento de uma bancada de testes para a verificação experimental dos diversos fenômenos de onda possíveis de serem simulados no laboratório de dinâmica de fluidos computacional. A bancada consiste de uma plataforma móvel, controlada por um motor elétrico, que, por sua vez, é comandado por um computador para seguir trajetórias de deslocamento previamente calculadas que gerem as ondas desejadas. No experimento são utilizadas interfaces eletrônicas que permitiram um controle de tipo digital implementado com auxilio do programa de desenvolvimento Labview. Os perfis de ondas são medidos através de fotografias e filmes, para posterior comparação com os perfis de ondas computacionais. Os resultados mostrados verificam uma boa coincidência entre ambos os perfis. Palavras-chave: Dinâmica dos fluidos computacionais, escoamento em superfície livre, controle digital.

ABSTRACT

Energetic waves may be formed when liquid storage tanks are shaked under certain

frequencies. These waves may cause structural damage or loss of maneuverability in vehicles.

In particular wave sloshing occurs when an energetic wave hits the walls of a tank. This

phenomenon has been studied numerically by some authors via the finite volume method

though its validation still relies on experiments. The present work aims to develop a shaking

table facility for wave sloshing analysis and validation of the numerical simulations. The

shaking table consists of a horizontal platform connected to an electrical motor, which is

controlled by a computer, leading to certain displacements and generating the desired free

surface waves at the reservoir. Electronic interfaces are used for digital control of the table

and implemented under the Labview environment. Free surface profiles are registered in the

form of photos and films for comparison with the simulated waves. Preliminary results show a

good agreement between the free surface profiles.

Keywords: Computational fluid dynamics, free surface flow, digital control.

DEDICATÓRIA

Este trabalho é dedicado a minha família: Tatiana Cristina Pocuve Saucedo e Lucas

Saucedo Ferreira; e aos meus pais Wilmar Dias Ferreira e Ana Maria Rangel.

Dedico este trabalho também as muitas pessoas que contribuíram para o alcance desse

objetivo, ajudando a gerar com sucesso a minha formação como mestre em engenharia

mecânica.

Espero que este trabalho possa contribuir para o desenvolvimento tecnológico no

Brasil, ajudando o país a se tornar um dos grandes líderes desenvolvedores de tecnologia no

mundo.

Espero que este trabalho possa ser para eu uma passagem para outro estágio da vida

profissional e pessoal.

AGRADECIMENTOS

Aos meus orientadores e professores Roger Matsumoto Moreira e José Andrés

Santisteban Larrea que me ajudaram no desenvolvimento deste trabalho, com empenho e

dedicação, transferindo com muito esforço o máximo de conhecimento possível para mim,

além de me dar todo o apoio e a oportunidade de poder estar fazendo este curso na

Universidade Federal Fluminense.

Ao programa de pós-graduação em engenharia mecânica da UFF pela oportunidade e

apoio institucional a mim dado no desenvolvimento dos meus estudos.

A UFF pela oportunidade que foi dada.

Aos alunos Elkin, Marcelo e Paolo por todo o suporte e toda a cooperação nas

simulações que envolveram o uso do software ANSYS CFX® e nos experimentos realizados

por este trabalho.

A todos os funcionários, professores e colegas da UFF.

Agradecimento em especial aos professores membros da banca e colaboradores:

Antônio Lopes Gama; Gustavo César Rachid Bodstein; Julio Tomás Aquije Chacaltana;

Domingos de Faria Brito David.

A todos os meus amigos que direta e indiretamente me apoiaram no desenvolvimento

desse trabalho.

EPÍGRAFE

Na vida existem apenas dois tipos de pessoas:

as que fazem as coisas acontecerem,

e as que falam sobre os que os outros fizeram.

Escolha em que lado você está e seja um vencedor.

Leonardo Ferreira, 2010

LISTA DE ILUSTRAÇÕES

Fig. 2.1 Um esboço do domínio fluido D, f. 17. Fig. 2.2 Histórico da pressão medida em uma parede vertical durante o impacto de uma

onda. (Peregrine 2004), f. 22. Fig. 2.3 Histórico da pressão em múltiplos locais para tanques (a) 30% e (b) 70%

cheios, f. 25. Fig. 2.4 Esquemas de interpolação temporal de uma propriedade genérica ψ, f. 34.

Fig. 3.1 Diagrama de blocos de um sistema de controle (a) em malha aberta, (b) em

malha fechada, f. 39.

Fig. 3.2 Diagrama de blocos de um sistema de controle digital em malha fechada, f. 40.

Fig. 3.3 Acionamento do motor em ponte completa, f. 42.

Fig. 3.4 Tensão no motor x tempo com uma estratégia PWM, f. 44. Fig. 3.5 Resposta ao degrau unitário de um sistema de primeira ordem, f. 46. Fig. 3.6 Diagrama de blocos de um sistema de controle em malha fechada, f. 47. Fig. 3.7 Lugar das raízes do sistema de controle com ganho proporcional, f. 48. Fig. 3.8 Lugar das raízes do sistema com controlador proporcional derivativo, f. 48. Fig. 3.9 Lugar das raízes do sistema com controlador proporcional derivativo e filtro, f.

49. Fig. 3.10 Lugar das raízes do sistema com controlador proporcional integral, f. 50. Fig. 3.11 Lugar das raízes do sistema com controlador proporcional integral derivativo,

f. 51. Fig. 4.1 Bancada de testes, f. 53. Fig. 4.2 Mesa móvel e reservatório, f. 54. Fig. 4.3 Ambiente de desenvolvimento, f. 56. Fig. 4.4 Interface gráfica do programa, f. 56. Fig. 4.5 Interface gráfica do programa, f. 57.

Fig. 4.6 Ambiente de desenvolvimento, f. 57. Fig. 4.7 Ambiente de desenvolvimento (cont.), f. 58. Fig. 4.8 Referência de deslocamento da mesa para o experimento N3, f. 59. Fig. 4.9 Lugar das raízes e pólos em malha fechada, f. 60. Fig. 4.10 Resposta experimental para o caso A, f. 61. Fig. 4.11 Resposta experimental para o caso B, f. 61. Fig. 4.12 Resposta experimental para o caso C, f. 62. Fig. 4.13 Resposta do controlador proporcional com ganho kp = 0,1 e referência de

posição com freqüência de 1,1 Hz, f. 62. Fig. 4.14 Resposta do controlador proporcional com ganho kp = 0,1 e referência de

posição com freqüência de 1,0 Hz, f. 63. Fig. 4.15 Três instantes em seqüência do movimento horizontal H4 do reservatório, f. 64. Fig. 4.16 Três instantes em seqüência do movimento horizontal H5 do reservatório, f. 65. Fig. 4.17 Pressão x tempo na parede do reservatório para o teste H4, f. 66. Fig. 4.18 Pressão x tempo na parede do reservatório para o teste H5, f. 66. Fig. 4.19 Pressão x tempo na parede do reservatório para o teste manual, f. 67. Fig. 5.1 Discretização empregada nas simulações com o software ANSYS CFX, f. 69. Fig. 5.2 Movimento harmônico lateral da mesa empregado no experimento H10 de

Bredmose et al. (2003). O deslocamento horizontal dipx é dados em metros e o tempo t em segundos, f. 69.

Fig. 5.3 Comparação entre os resultados numéricos (obtidos a partir do código

comercial ANSYS CFX) e experimentais (experimento H10 de Bredmose et al. 2003), f. 71/72.

Fig. 5.4 (a) Comparação entre as pressões medidas (- - - - -) e calculadas (_______) por

Bredmose et al. (2003) para o experimento H10. (b) Distribuição de pressões obtida através do código comercial ANSYS CFX, f. 73.

Fig. 5.5 Referência de deslocamento da mesa para o experimento N1, f. 74. Fig. 5.6 Três instantes em seqüência do movimento horizontal do reservatório, f. 75. Fig. 5.7 Interface gráfica do experimento N4, f. 76.

Fig. 5.8 Três instantes em seqüência do movimento do reservatório com o experimento N4, f. 77.

Fig. 5.9 Malha hexaédrica utilizada no experimento N1, f. 78. Fig. 5.10 À esquerda, resultados numéricos obtidos com o ANSYS CFX. À direita,

seqüência de fotos extraída do experimento N1, f. 79. Fig. 5.11 À esquerda, resultados numéricos obtidos com o ANSYS CFX. À direita,

seqüência de fotos extraída do experimento N4, f. 80/81. Fig. 5.12 Experimento N4: variação da pressão (a) medida e (b) calculada via CFX, f.

82.

SUMÁRIO

1 INTRODUÇÃO 14

1.1 MOTIVAÇÃO 14 1.2 OBJETIVO 15 1.3 SUMÁRIO 15

2 GERAÇÃO E PROPAGAÇÃO DE ONDAS EM UM RESERVATÓRIO 17

2.1 ESCOPO DO ESTUDO 17 2.2 PROBLEMA DE VALOR INICIAL 18 2.3 REVISÃO BIBLIOGRÁFICA 20 2.4 MÉTODO DE VOLUMES FINITOS 32 2.5 ANSYS CFX 35

3 MODELAGEM E CONTROLE 39

3.1 INTRODUÇÃO 39 3.2 CONTROLADORES 40 3.2.1 CONTROLADOR PROPORCIONAL (P) 41 3.2.2 CONTROLADOR INTEGRAL (I) 41 3.2.3 CONTROLADOR PROPORCIONAL-INTEGRAL (PI) 42 3.2.4 CONTROLADOR PROPORCIONAL-DERIVATIVO (PD) 42 3.2.5 CONTROLADOR PROPORCIONAL-INTEGRAL-DERIVATIVO (PID) 42 3.3 ACIONAMENTO COM ESTRATÉGIA DE CONTROLE PWM 43 3.4 MODELO DA PLANTA 45 3.5 PROJETO DO CONTROLADOR 47 3.5.1 CONTROLADOR PROPORCIONAL 47 3.5.2 CONTROLADOR PROPORCIONAL DERIVATIVO 48 3.5.3 CONTROLADOR PROPORCIONAL INTEGRAL 49 3.5.4 CONTROLADOR PROPORCIONAL INTEGRAL DERIVATIVO 50 3.5.5 DISCRETIZAÇÃO DO CONTROLADOR CONTÍNUO 52

4 DESCRIÇÃO DA BANCADA DE TESTES 53

4.1 COMPONENTES 53 4.2 PROGRAMAS DE CONTROLE DA MESA 55 4.3 EXPERIMENTO N3: TESTE EM MALHA FECHADA 59

5 RESULTADOS 68

5.1 VALIDAÇÃO VIA EXPERIMENTO H10 (BREDMOSE ET AL. 2003) 68 5.2 EXPERIMENTO N1 74 5.3 EXPERIMENTO N4 75 5.4 SIMULAÇÃO NUMÉRICA E VALIDAÇÃO VIA EXPERIMENTO N1 78 5.5 SIMULAÇÃO NUMÉRICA E VALIDAÇÃO VIA EXPERIMENTO N4 79

6 CONCLUSÕES 83

7 REFERÊNCIAS BIBLIOGRÁFICAS 85

8 APÊNDICE 94

PROJETO DO CONTROLADOR PI 94 DISCRETIZAÇÃO DO CONTROLADOR CONTÍNUO 96

14

1 Introdução 1.1 Motivação

Os desenvolvimentos da eletrônica e da capacidade de processamento dos

computadores permitiram, nas últimas décadas, que simulações computacionais através de

métodos numéricos pudessem representar experimentos que demorariam meses ou anos para

serem desenvolvidos. Em particular, muitos processos relacionados a fenômenos de transporte

podem ser simulados hoje através de métodos computacionais, reduzindo de uma forma

significativa o tempo de desenvolvimento de projetos experimentais. Contudo, protótipos que

reproduzam em pequena escala esses fenômenos são imprescindíveis para a validação das

simulações numéricas. Felizmente, o desenvolvimento da eletrônica de potência evoluiu a tal

ponto que é possível reproduzir uma grande gama de fenômenos na área de controle de

processos.

No âmbito da mecânica dos fluidos, um fenômeno de particular interesse é o

estudo da geração e propagação de ondas em meios confinados. A agitação de um reservatório

parcialmente cheio de líquido pode gerar ondas gravitacionais extremamente energéticas. Este

fenômeno, conhecido por sloshing, é um problema de autovalor clássico na Mecânica dos

Fluidos. Grandes cientistas, como Poisson, Rayleigh e Kirchhoff, já realizaram investigações

sobre este tema (para uma revisão completa, veja Lamb 1932) e, atualmente, vários livros

tratam especificamente sobre esse assunto (Ibrahim 2005, Faltinsen & Timokha 2009).

Tanques de combustível em caminhões, aviões, espaçonaves e navios podem, sob

certas circunstâncias, estarem sujeitos a cargas harmônicas de alta amplitude, podendo

alcançar modos de ressonância. Tais ondas, quando refletidas nas paredes do compartimento,

impõem forças hidrodinâmicas que podem gerar um efeito desestabilizante. Esta interação

entre o fluido e a estrutura pode causar a ruptura do container, instabilidade e perda de

manobrabilidade em veículos que transportam substâncias líquidas. Navios petroleiros e

plataformas do tipo FPSO (Floating Production, Storage and Offloading Unit) e FSO

(Floating Storage and Offloading Unit) são exemplos de veículos que transportam tanques de

carga de grande porte e que naturalmente estão sujeitos aos movimentos de jogo,

afundamento, arfagem, deriva, avanço e guinada.

É de suma importância conhecer o mecanismo de formação dessas ondas,

diagnosticando os principais modos de excitação do sistema que poderiam acarretar em

impactos violentos de líquido nas estruturas do container e comprometer a estabilidade do

15

veículo transportador. Neste sentido, pesquisadores da área têm atuado na identificação das

freqüências excitadoras responsáveis pelo sloshing e na modelagem da propagação dessas

ondas em meios confinados (Bredmose et al. 2003, Moreira et al. 2006, Mendes 2007, Kim

2007, Nicolato 2009, Nicolato & Moreira 2009). A magnitude das pressões exercidas por

essas ondas sobre as paredes do reservatório também são de interesse do ponto de vista

científico (Peregrine 2004), além de servir como base para o dimensionamento de vasos e

tanques de armazenamento de óleo e GNL (gás natural liquefeito). Atualmente, projetistas

têm utilizado códigos de fluidodinâmica computacional e experimentos na tentativa de

entender melhor o processo de formação dessas ondas e os impactos causados nas estruturas

do reservatório.

1.2 Objetivo

O presente trabalho tem como objetivo desenvolver uma bancada de testes para a

verificação experimental da formação de ondas em um reservatório prismático. A bancada

consiste de uma plataforma móvel controlada por um motor elétrico que, por sua vez, é

comandada por um computador para seguir determinadas trajetórias de deslocamento que

gerem as ondas desejadas no reservatório. No experimento, são utilizadas interfaces

eletrônicas que permitem um controle do tipo digital, implementado com auxílio do programa

de desenvolvimento Labview. Os perfis da superfície livre são registrados na forma de fotos e

filmes, para posterior comparação com as ondas simuladas numericamente. As pressões

exercidas pelas ondas nas paredes do reservatório são medidas com a ajuda de um sensor e

posteriormente são processadas.

O movimento não-estacionário do escoamento com superfície livre é simulado

numericamente através do código de fluidodinâmica computacional ANSYS CFX. O

problema de valor de contorno é resolvido via método de volumes finitos com um modelo

multifásico homogêneo para a superfície livre. Os resultados numéricos são comparados aos

experimentais em termos dos perfis de ondas gerados e das pressões exercidas pelas ondas

sobre as paredes do reservatório.

1.3 Sumário

A estrutura da dissertação é a seguinte. O capítulo 2 define o escopo do estudo, as

hipóteses assumidas na modelagem matemática, suas equações governantes e condições de

contorno. Uma revisão bibliográfica das principais técnicas numéricas empregadas na solução

desse problema de valor de contorno e algumas validações experimentais é mostrada, bem

16

como o método numérico empregado na simulação do tanque de ondas. O capítulo 3

apresenta os fundamentos da teoria de controle e alguns dos tipos de controladores mais

utilizados no mundo industrial. O capítulo finaliza com a apresentação do modelo matemático

da planta em estudo e o projeto do controlador utilizado. No capítulo 4 são apresentados os

detalhes da bancada experimental, incluindo os circuitos eletrônicos auxiliares e os programas

computacionais desenvolvidos. No capítulo 5 são analisados os resultados experimentais

encontrados, comparando-os com as simulações numéricas correspondentes. São monitorados

os perfis da superfície livre e as pressões exercidas em uma das paredes do reservatório. Por

fim, o capítulo 6 apresenta as principais conclusões e sugestões para trabalhos futuros.

17

2 Geração e Propagação de Ondas em um Reservatório

2.1 Escopo do estudo

Este trabalho tem por objetivo modelar numérica e experimentalmente o

escoamento com superfície livre que é formado devido à perturbação de um reservatório

parcialmente cheio de líquido sujeito à ação do campo gravitacional gr . A figura 2.1 ilustra

um plano secante x-y do domínio fluido tridimensional, com suas coordenadas

correspondentes e parâmetros. Por limitações do aparato experimental, apenas movimentos

laterais na direção x são considerados em nosso estudo. Assim, toda a informação necessária

para a descrição do escoamento está contida à priori no plano x-y, contendo o fluido e sua

superfície livre, sendo delimitados por duas paredes verticais e um fundo horizontal, todos

impermeáveis.

Figura 2.1. Um esboço do domínio fluido D.

O deslocamento da superfície livre é medido a partir de uma distância da condição

de repouso absoluto. Coordenadas cartesianas são definidas com o eixo x na superfície livre

em repouso, com o eixo y apontando para cima (de tal forma que o fluido ocupe o plano y ≤ 0

quando em repouso) e com o eixo z perpendicular ao plano secante x-y. Dois reservatórios em

vidro com dimensões internas de 393 x 194 x 265 mm³ (comprimento x largura x altura) e

800 x 100 x 400 mm³ são empregados nos experimentos (para maiores detalhes, vide capítulo

5). Em nosso modelo numérico, o escoamento é considerado incompressível e o fluido

viscoso, enquanto os efeitos da tensão superficial são desprezados. Por simplicidade,

consideraremos a água como fluido de trabalho. A motivação para tal estudo encontra-se na

18

seção 1.1. O problema de valor inicial a ser resolvido numericamente é apresentado na

próxima seção.

2.2 Problema de valor inicial

As equações de conservação da massa e da quantidade de movimento devem ser

satisfeitas no domínio fluido D. Considerando um volume V de fluido, com massa específica

ρ, a massa dentro do volume de controle é dada por ∫V

dVρ . A taxa com que a massa flui

através da superfície de controle pode ser descrita como ∫ ⋅A

dAnu rrρ , onde ur=(u,v,w) é o

campo de velocidades e nr é o vetor normal à superfície de controle. dV e dA são elementos

infinitesimais do volume confinado e da área da superfície de controle.

Para um modelo bifásico, a fração volumétrica ϕ de cada fluido (água, ϕ1, e ar, ϕ2)

em cada volume de controle define as propriedades do fluido. O ar é considerado um gás ideal

e supõe-se um escoamento incompressível para a água. O balanço de massa na forma

diferencial é então descrito por:

( ) ( ) 0=⋅∇+∂

∂ ut

rφρφρ . (2.1)

A conservação da quantidade de movimento é regida pela equação de Navier-

Stokes:

( ) fTDt

uD rrrrρφφρφ

+⎟⎠⎞⎜

⎝⎛⋅∇= , (2.2)

onde Trr

é o tensor das tensões (representando a soma das tensões normais, viscosas e de

expansão volumétrica) e fr

ρ as forças de corpo atuantes no volume V. Consideraremos que

apenas a aceleração gravitacional atue externamente ao volume de controle, tal que

( ))0,,0 gf −=r

.

Para um fluido newtoniano com viscosidade μ, viscosidade de expansão k e

submetido a uma pressão dinâmica p, o tensor Trr

pode ser escrito em notação indicial como:

19

( ) ⎟⎟⎠

⎞⎜⎜⎝

⎛

∂∂

+∂∂

+∂∂

−+−=i

j

j

iij

i

iijij x

uxu

xukpT μδμδ 3

2 . (2.3)

Para um modelo homogêneo bifásico:

∑=

=2

1lllρφρ , ∑

=

=2

1lllμφμ . (2.4)

Como condições de contorno do problema, consideram-se o fundo horizontal e as

paredes verticais rígidas e impermeáveis. A condição de não-deslizamento é também aplicada

em todas as superfícies do reservatório.

),0,0(),( == RRR vuur (2.5)

onde Rur é a velocidade relativa do fluido em relação à parede do reservatório.

Como condição inicial do problema, supõe-se o estado de repouso absoluto para a

superfície livre para t=0, imediatamente seguido de uma perturbação senoidal do reservatório

na direção x:

x(t) = x0 sen (ωt), (2.6)

onde x0 é a amplitude do movimento e ω a frequência de oscilação do reservatório.

Nos modelos de resolução de problemas de transferência de calor e massa, a

turbulência exerce normalmente um papel importante e pode ser tratada matematicamente

através da proposição de um modelo de fechamento das equações médias espaciais e

temporais de Navier-Stokes, também denominado de RANS (Reynolds Averaged Navier-

Stokes Equations). Entre os mais difundidos métodos nessa classe, e que usam duas equações

diferencias parciais adicionais, estão os modelos k-ε e o SST (Shear Stress Transport).

O modelo SST calcula o transporte da tensão cisalhante turbulenta e garante boa

precisão na estimativa da quantidade de separação do escoamento submetido a gradientes de

20

pressão adversos. Largamente utilizado em representação de escoamentos externos. Esse

modelo foi empregado nas simulações realizadas por este trabalho.

2.3 Revisão bibliográfica

Na Mecânica dos Fluidos, o fenômeno sloshing representa um problema de

autovalor clássico, tendo sido objeto de estudo de cientistas ilustres como Poisson, Rayleigh e

Kirchhoff. Uma referência padrão sobre o assunto é Lamb (1932, §258), que relata a escassez

de soluções para a superfície livre supondo um reservatório com profundidade variável. De

fato, Lamb analisou o movimento do fluido em um canal triangular, cuja seção transversal

consistia de dois segmentos de reta inclinados a 45o em relação à vertical, que permanece

como uma das poucas soluções analíticas encontradas até hoje.

Em meados do século XX, houve mais uma vez o interesse da comunidade

acadêmica pelo desenvolvimento de estudos que envolviam problemas associados a tanques

de combustível em foguetes e aviões. Em um tanque bidimensional com comprimento πL e

profundidade hL, Moiseyev (1958) e Chester (1968) identificaram que as oscilações na

superfície livre possuem um espectro discreto de freqüências ( )[ ] 21/tanh Lnhng=ω , onde g

é a aceleração da gravidade e n é o número de ondas.

Ockendon & Ockendon (1973) constataram que pequenas oscilações verticais e

horizontais de um container parcialmente cheio de líquido podem gerar uma movimentação

do fluido em uma escala muito maior do que a amplitude de excitação do sistema. De fato,

por intermédio de métodos assintóticos, Ockendon & Ockendon reproduziram os resultados

analíticos de Moiseyev (1958) e Chester (1968), além de unificar sua análise para águas rasas

(h/L<1). Moiseev (1964) e Moiseev & Petrov (1968) fornecem uma revisão extensa do

fenômeno sloshing à luz da teoria linear, com soluções via teoria da perturbação, não tratando,

no entanto, de oscilações não-lineares.

Ainda no campo analítico, Davis (1965) obteve, através de métodos assintóticos,

os autovalores de oscilações bidimensionais em canais de seção transversal arbitrária. Em dois

artigos seguintes, Davis (1974) estendeu a metodologia anterior a geometrias semi-circulares

e Packham (1980) resolveu o caso de um canal de seção triangular, supondo uma inclinação

de 30o em relação à horizontal.

Outras referências clássicas no estudo de ondas ressonantes em ambientes

confinados incluem Ursell (1952), que analisou a formação de ondas laterais em um container

prismático, Benjamin & Ursell (1954), que considerou oscilações verticais de um tanque

cilíndrico de seção transversal arbitrária, e Faltinsen (1974), que encontrou soluções analíticas

21

para a superfície livre imprimindo pequenas movimentações de jogo e deriva ao reservatório

de seção retangular.

Em um trabalho analítico recente, Chapman & Porter (2005) resolveu o problema

de sloshing bidimensional em um tanque de paredes laterais verticais e fundo arbitrário. Para

determinação da primeira e segunda ordem do potencial de velocidades na superfície livre,

empregaram-se as funções de Green para um leito plano, em conjunto com as equações de

Cauchy-Riemann para um fundo variável.

Com o aparecimento de computadores com maior capacidade de processamento e

os crescentes avanços na área de análise numérica abriram-se novos horizontes no estudo de

escoamentos com superfície livre, onde efeitos não-lineares desempenham um papel

importante. De fato, o desenvolvimento de novas técnicas de fluidodinâmica computacional

tornou possível a investigação de determinados fenômenos através de meios que não sejam

exclusivamente experimentais.

Longuet-Higgins & Cokelet (1976) foram os precursores na aplicação do método

da integral de contorno na modelagem de ondas não-lineares, seguidos por uma série de

outros autores que estenderam tal técnica para computar outras configurações de escoamento

com superfície livre (ver, por exemplo, Vinje & Brevig (1981), Baker et al. (1982), Dold &

Peregrine 1986). Atualmente, o método da integral de contorno é aplicável às mais diferentes

geometrias, em domínios confinados ou semi-infinitos, fornecendo soluções para ondas

gravitacionais e capilares em regimes estacionários ou transientes.

Mercer & Roberts (1992, 1994) foram os primeiros a modelarem de forma precisa

ondas estacionárias íngremes bidimensionais. Eles verificaram que ondas estacionárias muito

íngremes em águas profundas são afetadas por instabilidades numéricas. Outros estudos

envolvendo a evolução de ondas completamente não-lineares incluem McIver & Peregrine

(1981), Anderson et al. (1990), Topliss (1994), Longuet-Higgins (2000, 2001a) e Longuet-

Higgins & Dommermuth (2001a,b).

Wu & Eatock Taylor (1994) aplicaram o método dos elementos finitos ao

problema de sloshing em um tanque bidimensional, obtendo soluções para o movimento

fluido pela aproximação de Galerkin. Em um trabalho seguinte, Wu & Eatock Taylor (1998)

estenderam seu modelo para o problema 3D, onde foram observadas, além das ondas

estacionárias, ondas progressivas e ressaltos hidráulicos

Em um trabalho extenso, Bredmose et al. (2003) analisaram experimentalmente e

numericamente a geração de ondas não-lineares em um reservatório de dimensões 1480 x 400

x 750 mm3 (comprimento x largura x altura). Oscilações horizontais são impostas com o

22

objetivo de se investigar o impacto de ondas nas paredes laterais do tanque. Neste caso, um

modelo estendido das equações de Boussinesq é utilizado, seguindo Wei et al. (1995), e uma

boa concordância entre os resultados numéricos e experimentais é encontrada, excetuando-se

na região próxima à parede vertical, onde o modelo superestima o run up e não reproduz o

jato descendente. Tais diferenças são creditadas à falha do modelo empregado, que considera

como condição de contorno um leito com batimetria suave.

Impactos violentos de ondas gravitacionais podem resultar em danos ou colapsos

estruturais e sempre representaram uma preocupação para os projetistas. Em uma larga escala,

a distribuição de pressão quando uma onda gerada por uma tempestade encontra um quebra-

mar tem sido objeto de discussão de vários artigos (Chan 1994, Zhang et al. 1996). Estudos

experimentais (Hattori et al. 1994) e teóricos (Cooker & Peregrine 1992) têm focalizado nas

pressões exercidas pelo impacto de ondas em paredes verticais. A fim de ilustrar este

fenômeno a figura 2.2 mostra o perfil típico da pressão medida em uma parede vertical

durante o impacto de uma onda (Peregrine 2004). Pode-se notar a formação de dois picos, o

primeiro predominantemente devido à inércia do fluido, visto que a escala de tempo (da

ordem do milissegundo) é muito pequena para que a gravidade exerça qualquer influência, e o

segundo devido à pressão hidrostática gerada pelo movimento descendente da onda após o

máximo run up.

Figura 2.2. Histórico da pressão medida em uma parede vertical durante o impacto de

uma onda. (Peregrine 2004)

Muitos modelos matemáticos têm sido empregados na previsão do problema de

sloshing. Por exemplo, a teoria linear para escoamentos potenciais provê meios simplificados

de avaliar as cargas induzidas pela geração de ondas e oferece critérios físicos do mecanismo

23

do fenômeno, mas é utilizado de forma restrita. Para casos envolvendo tanques com níveis de

preenchimento baixos, onde a razão da profundidade do líquido em relação às dimensões do

tanque é menor que 0,20, para geometrias complexas ou quando obstáculos internos são

usados, a teoria potencial é falível.

Normalmente, para escoamentos rasos em tanques as equações de Navier-Stokes

ou suas formas simplificadas, como as equações de Euler onde a viscosidade da fase líquida é

negligenciada, são usadas para estimar o fenômeno em estudo. O cálculo do perfil da

superfície livre em um tanque é um elemento chave na aproximação precisa das cargas

geradas pelo sloshing. O método do marcador e célula (marker and cell method, MAC),

método do volume de fluido (volume of fluid method, VOF), hidrodinâmica de partículas

suavizadas (smooth particle hydrodynamics, SPH) e o método de ajuste de nível (level set

method, LSM) ou a combinação do LSM com o VOF têm sido aplicados com sucesso para o

cálculo da superfície livre.

Ainda no campo numérico / experimental, Kim (2007) reporta as observações de

escoamentos de sloshing fortemente não-lineares em navios cargueiros e seus efeitos

acoplados com o movimento da embarcação. Dois diferentes modelos numéricos são

utilizados:

• O primeiro esquema se baseia na discretização do domínio em uma malha

finita, adotando o esquema colocalizado para o cálculo do perfil de

velocidades, a partir das equações de Euler e da continuidade – Método das

diferenças finitas.

• O segundo esquema se baseia no conceito da partícula onde as equações do

domínio são resolvidas individualmente para cada partícula, levando-se em

consideração uma função peso que reflita a influência das partículas vizinhas –

Método SPH.

Kim (2007) reporta diversos aspectos físicos do sloshing observado nos

experimentos realizados. Fenômenos do escoamento violento como quebra de ondas,

esguichos (splashes), e impactos são comumente observados no sloshing em cargueiros. O

impacto de pressão decorrente da geração de ondas é afetado por estes fenômenos, tornando a

observação e o entendimento destes fenômenos físicos importantes para o desenvolvimento de

modelos numéricos adequados. Para a observação dos fenômenos listados a seguir, um tanque

menor que o utilizado pelos experimentos de Bredmose et al. (2003)é submetido a oscilações

laterais (sway).

24

• Formação de bolsões de ar

Bolsões de ar podem ser formados nas extremidades do teto do tanque ou nas

laterais durante um impacto, e desempenham papel de amortecimento dos

picos de pressão decorrentes do impacto. Os resultados dos experimentos

revelam que quando um tanque tem níveis de preenchimento elevado, a

formação de bolsões de ar não se dá apenas nas bordas superiores do tanque,

mas também na parte central do teto do tanque.

• Formação de bolhas de gás

As bolhas são formadas devido à entrada de gás no fluido durante movimentos

bruscos da superfície livre, por exemplo, durante a quebra da onda ou durante

o impacto. Um aglomerado de bolhas reduz a massa específica do fluido, o que

explica a redução do impacto de pressão.

• Quebra de ondas e propulsão de esguichos

Em geral, a não-linearidade do sloshing se torna dominante para tanques pouco

preenchidos. Para escoamentos rasos, a formação de um ressalto hidráulico ou

quebra de onda é frequente. Um típico processo de impacto em condições rasas

é a formação da superfície livre fortemente não-linear, seguida pelo impacto

nas paredes laterais e run-up nas paredes, com possíveis splashes. Por outro

lado, a quebra de ondas em tanques mais cheios, não ocorre com freqüência,

mas uma grande quantidade de esguichos é formada após impactos no topo do

tanque.

• Hidroelasticidade

A figura 2.3 mostra dois exemplos de um histórico de pressão medido em

múltiplos locais nas paredes do tanque num mesmo momento. Após a

ocorrência do impacto nas paredes laterais (a) e no topo (b) do tanque, os sinais

de pressão mostram oscilações de grande freqüência, o que pode ser explicado

pela hidroelasticidade das estruturas do tanque. Considerando o período e a

magnitude da oscilação, é evidente a vibração da parede lateral após o impacto

no topo do tanque.

25

Figura 2.3. Histórico da pressão em múltiplos locais para tanques (a) 30% e (b) 70%

cheios.

Outras referências de trabalhos onde a geração de ondas é modelada

numericamente incluem Kyoung et al. (2005) que simulam numericamente a carga do

impacto de onda devido a um violento sloshing através do método dos elementos finitos

aplicados às condições não-lineares de superfície livre. São desenvolvidas as equações

governantes do domínio fluido aplicadas à formulação variacional. A computação é

processada com o objetivo de prever o impacto da onda em um pilar localizado no meio de

um tanque sob oscilações laterais harmônicas. É observado numericamente que o impacto da

onda aumenta na medida em que a profundidade da água aumenta, para profundidades baixas

relativas ao comprimento do tanque. Cox et al. (2003) estudaram os efeitos da mudança

gradual no comprimento do tanque em ondas estacionárias não-lineares e em ondas

ressonantes em tanques rasos com água. A análise envolve o desenvolvimento de modelos

numéricos a partir das equações de Boussinesq.

Ma & Yan (2005) apresentam um novo modelo numérico para a simulação do

sloshing baseado na teoria não-linear, denominado quasi arbitrary Lagrangian-Eulerian finite

element method (QUALE-FEM). A principal diferença do método apresentado para o método

de elementos finitos tradicional se dá na geração da malha computacional que é feita apenas

uma vez no início da simulação e é adaptada nos passos de tempo subsequentes para

conformar com o movimento da superfície livre e com as estruturas. O modelo é validado por

ondas geradas em um tanque de água com resultados satisfatórios quando comparados com

soluções analíticas, dados experimentais e outros modelos numéricos.

Um modelo de elemento de contorno é apresentado para a análise do movimento

acoplado da vibração de estruturas com sloshing de pequenas amplitudes em um domínio

26

bidimensional por Yuelong & Xing (1999) em um trabalho onde as equações de Navier-

Stokes são linearizadas e transformadas em equações de integral de contorno.

Rabier & Medale (2003) apresentam um modelo numérico para a simulação de

escoamentos viscosos e incompressíveis com superfície livre, que se baseia na resolução das

equações governantes através do método dos elementos finitos. Uma malha móvel é

introduzida para lidar com a presença da superfície livre e a validação é realizada com

problemas de sloshing bidimensional com pequenas deformações.

Chen et al. (2008) investigam numericamente a dinâmica do impacto de ondas em

um tanque oscilando horizontalmente de forma harmônica com diferentes níveis de

preenchimento e freqüências de excitação, baseando-se no método do ajuste de nível (LSM).

Investigam também a dependência da solução numérica na resolução da malha, tamanho do

passo temporal e a espessura da interface. Para os casos de nível de preenchimento elevado o

esquema numérico de predição do impacto de pressão provê níveis de precisão aceitáveis

comparando-se com dados experimentais. Para níveis de preenchimento baixo a mesma

concordância foi observada, porém os esquemas numéricos demandaram uma discretização

elevada da malha espacial e temporal, conseqüentemente um aumento do recurso

computacional. Testes também são conduzidos para tanques retangulares com chanfros no

topo (45° e 60°) sujeitos a oscilações angulares, onde se observou menores valores máximos

de pressão para os chanfros de 45°.

Godderidge et al. (2008) publicaram um trabalho onde foram exploradas as

diferenças entre os modelos multifásicos homogêneo e não-homogêneo. Para isso, um

escoamento em um reservatório com água oscilando próximo às freqüências de ressonância

foi usado.

A condição multifásica homogênea utiliza das frações volumétricas de cada fase

do domínio fluido para determinar suas propriedades em cada célula da malha espacial. É

também conhecido como método do volume de fluido. Outra abordagem mais geral, porém

computacionalmente mais cara, é o modelo multifásico não-homogêneo, onde a solução dos

campos de velocidade é feita separadamente para cada fluido, combinada nas interfaces, com

modelos de transferência de massa e momento. Devido ao baixo custo computacional e boa

estabilidade numérica, o modelo homogêneo é preferencialmente utilizado nos estudos de

sloshing. Para situações onde há grande interação entre as fases, o segundo modelo se mostra

mais adequado.

Godderidge et al. concluem, a partir de análises dimensionais dos resultados

numéricos, que o modelo não-homogêneo é mais apropriado para situações de sloshing

27

violento. Concluem que a velocidade relativa entre as fases, que é ignorada pelo modelo

homogêneo, pode ser significante se comparada ao campo global de velocidade. O histórico

de pressão do modelo não-homogêneo verificou bem os resultados experimentais, enquanto

que o histórico de pressão do modelo homogêneo, para as simulações realizadas, apresenta

algumas deficiências, subestimando os dados experimentais. Por outro lado, o modelo não-

homogêneo apresentou-se de duas a três vezes mais caro computacionalmente.

Quando se compara escoamentos do sloshing em condições de temperatura e

pressão similares, em freqüências de oscilação próximas às de ressonância, o movimento da

superfície livre depende basicamente de três parâmetros que são o amortecimento devido à

viscosidade, a freqüência de oscilação associada à amplitude do movimento e a profundidade

do fluido. Royon-Lebeaud et al. (2006), com o objetivo de estudar a influência das

freqüências de oscilação em reservatórios cilíndricos, realizaram experimentos com líquidos

poucos viscosos, e com grande profundidade de forma que estas não se tornassem um

parâmetro de análise. Os limites da existência de diferentes regimes de ondas (denominadas

como ondas planares, espirais, sloshing caótico e quebra de ondas) foram determinadas como

função das freqüências motrizes em relação à menor freqüência natural, para uma variedade

de amplitudes.

Os estudos revelaram que perfis de onda planar em reservatórios cilíndricos são

observados para freqüências excitantes ligeiramente maiores que a menor freqüência natural

de oscilação, podendo bifurcar em ondas espirais (redemoinhos) para uma amplitude

limítrofe, que depende do valor da freqüência. Sloshing caótico e quebra de ondas ocorrem

quase periodicamente, iniciando-se com o crescimento da amplitude do movimento das ondas

planares seguido geralmente por quebra e sloshing irregular, então, após um amortecimento

do movimento irregular, a amplitude do movimento das ondas planares aumenta novamente.

Este comportamento é previsto para tanques cilíndricos, com grande profundidade de fluidos

poucos viscosos de acordo com os resultados experimentais levantados.

Cargueiros de GNL com tanques do tipo membrana são atualmente projetados

com um acréscimo significante na capacidade total de carga e nos tamanhos. Terminais de

carregamento e descarregamento não oferecem restrições para os níveis de enchimento de

GNL nos tanques. Estimativas precisas dos impactos do movimento da superfície livre interna

nas estruturas do tanque são necessárias. Daí o aparecimento recente de estudos do sloshing

aplicados ao armazenamento do gás natural liquefeito.

O instituto de pesquisas norueguês MARINTEK, parceira do CeSOS (Center of

Ships and Offshore Structures), desempenhou um grande número de testes para diferentes

28

geometrias de cargueiros e tanques do tipo membrana. O principal objetivo destes testes é

identificar propriedades estatísticas dos picos de pressão a respeito da magnitude, local e

duração. O tanque modelo é submetido a movimentos irregulares, simulando situações

marítimas reais.

A base para a relação de testes é identificar combinações críticas de níveis de

enchimento e parâmetros de excitação do tanque. A influência do sloshing no movimento dos

navios também é investigada nos experimentos realizados.

Chun et al. (2008) investigam experimentalmente a força dinâmica característica

dos sistemas de carga de GNL sobre cargas de impacto. Uma série de testes de impacto em

sistemas de isolamento de tanques tipo membrana de GNL foram realizados.

Graczyk & Moan (2008) investigam uma grande amostra de medidas de impacto

de pressão decorrente do sloshing de GNL, focado na magnitude dos impactos individuais e

suas configurações temporal e espacial. Os experimentos foram conduzidos em tanques com

níveis de enchimento altos (92.5%) e baixos (30%). A abordagem de pós-processamento de se

considerar um perfil triangular para o histórico de pressões é utilizada, e a aproximação para

um perfil trapezoidal é proposta.

Os experimentos foram baseados em modelos de tanques do tipo membrana de

cargueiros de GNL. O movimento aplicado considera condições de mar reais. Foi observado

que, enquanto impactos mais fracos decorrentes do sloshing podem produzir históricos

temporais complexos, o número de distorções locais (como picos secundários) decrescem

significativamente para os impactos de pressões máximas mais elevadas.

Comparando-se os perfis do histórico de pressão reais e aproximados (triangular),

foi observado que os tempos de ascensão não foram modelados adequadamente, daí a

proposta da abordagem do perfil trapezoidal. As características da magnitude de pressão e os

perfis temporais e espaciais dos impactos decorrentes do sloshing foram analisados

separadamente.

As pressões máximas registradas nos experimentos com nível de preenchimento

baixo foram maiores do que as registradas nos níveis de preenchimento alto. Quando todos os

históricos temporais de magnitude do impacto são considerados, a duração e a distribuição do

tempo de ascensão registram pequenos impactos com históricos mais complexos, o que não é

registrado para as amostras com nível de preenchimento alto.

Os impactos para altos níveis de preenchimento são mais localizados

espacialmente, além disso, a razão entre pressão média espacial e a pressão máxima aumenta

significativamente quando os impactos com maiores pressões máximas são considerados. Os

29

impactos decorrentes do sloshing em tanques com níveis de preenchimento baixo geralmente

têm maior magnitude, durações maiores, menor razão entre o tempo de ascensão e duração e

maior extensão espacial. Todos estes fatores poderiam, a princípio, aumentar a resposta da

estrutura. A direção dominante do perfil de pressão é vertical para baixos níveis de

preenchimento e longitudinal para altos níveis de preenchimento.

Lee et al. (2006) desenvolveram um estudo de sensibilidade em determinados

parâmetros do sloshing em tanques de gás natural liquefeito utilizando um programa de

fluidodinâmica computacional baseado no método dos volumes finitos (FLOW-3D).

Primeiramente, uma breve análise dimensional é conduzida com o objetivo de identificar os

parâmetros governantes do escoamento, a qual indica dez parâmetros adimensionais que

governam os picos de pressão decorrentes do impacto de ondas. Dentre estes, quatro

parâmetros se referem a características geométricas do reservatório e um deles denota o

número de Froude (definido como a razão das forças de inércia e gravitacionais). Quando este

último é mantido em testes, a similaridade depende dos demais parâmetros:

• Número de Reynolds (razão entre as forças de inércia e forças viscosas)

• Razão entre a massa específica das fases vapor e líquida

• Pressão da fase vapor

• Compressibilidade das fases vapor e líquida

• Elasticidade da parede

O estudo foca nos efeitos da turbulência, viscosidade dos fluidos, razão entre

massas específicas, pressão da fase vapor e compressibilidade no impacto de pressão do

sloshing de GNL. Os efeitos da variação da viscosidade no impacto da pressão se mostraram

negligenciáveis de acordo com os resultados, o que levanta um ponto interessante a ser

analisado, visto que o estudo anterior, Royon-Lebeaud et al. (2006), releva justamente a

influência desta no fenômeno em questão.

Como pode ser observado nos estudos reportados, a geração de ondas em

reservatórios pode acarretar em forças hidrodinâmicas acentuadas nas paredes de um tanque,

o que pode ser um perigo para a integridade estrutural do equipamento. Com o objetivo de

minimizar estas cargas, para se prevenir falhas estruturais e se almejar um controle da posição

do fluido no tanque, extensos estudos teóricos e práticos têm sido elaborados tendo como

conseqüência o desenvolvimento de dispositivos que possam ser projetados para suprimir ou

mitigar a mobilidade do fluido.

30

Anderson et al. (2000) propõem a introdução de absorvedores de sloshing como

uma alternativa prática para os sistemas de controle de amortecimento. Adicionalmente, um

procedimento numérico de simulação é introduzido. Absorvedores de sloshing são osciladores

mecânicos cuja freqüência de ressonância é ajustada a uma freqüência crítica da estrutura a

ser controlada, e são capazes de contabilizar a energia dissipada pelo sloshing.

Geralmente, os efeitos decorrentes da geração de ondas são suprimidos de formas

passivas, através da introdução de subestruturas adicionais nos tanques. A introdução de

chicanas altera a freqüência de oscilação natural do fluido. Abramson (1969) reduziu os

parâmetros primários que influenciam o projeto dos dispositivos de supressão do sloshing

como:

• o projeto do veículo que transporta o tanque e sua trajetória;

• a amplitude do movimento do sloshing para diferentes níveis de

preenchimento;

• as características físicas do tanque, como geometria, deformação elástica e

isolamento;

• o nível de preenchimento e necessidade de drenos;

• o impacto que pode ser sustentado pelos dispositivos.

O conceito básico da proteção passiva é dissipar a energia do movimento através

da quebra do escoamento em diversas subcorrentes. O critério de projeto da NASA sugere a

subdivisão dos tanques longitudinalmente por paredes verticais ou a instalação de obstáculos

(chicanas) na direção da propagação da onda. A inserção de obstáculos no sentido da

propagação da onda é também um tópico apresentado em handbooks de construção offshore

para se prevenir o impacto decorrente do sloshing.

Gavrilyuk et al. (2005) estuda as soluções fundamentais do problema em um

tanque cilíndrico vertical com finas chicanas rígidas horizontais. São desenvolvidos métodos

analíticos para o problema linearizado de sloshing que aproxima as freqüências naturais e o

potencial de velocidades de forma precisa. O amortecimento decorrente da inserção das

chicanas é quantificado.

Cho et al. (2005) analisam numericamente as características de ressonância em

um tanque bidimensional com chicanas laterais (paralelas à superfície livre em seu estado

inicial). Um termo artificial de amortecimento é empregado na condição de contorno de

superfície livre para refletir os efeitos dissipativos no sloshing ressonante baseado no método

numérico de elementos finitos. Os resultados são validados com modelos analíticos para o

caso do tanque sem chicanas.

31

Os efeitos do nível de preenchimento do tanque, o número de chicanas, o

comprimento de abertura das chicanas e a localização destas foram parametricamente

analisados. A freqüência de ressonância e a altura máxima atingida pela superfície livre se

mostraram decair uniformemente com o aumento do número de chicanas, da altura de

instalação destas, da redução do comprimento de abertura destas e do nível de preenchimento

do líquido. O mesmo decréscimo com o aumento dos parâmetros citados ocorre para a pressão

hidrodinâmica e a força dinâmica resultante. Entre os quatro parâmetros citados, a influência

da altura da instalação das chicanas se mostrou mais relevante.

Um método multidisciplinar de projeto e otimização é apresentado por Craig &

Kingsley (2006) para dar suporte aos projetos de reservatórios parcialmente cheios. O código

comercial de CFD (Computational Fluid Dynamics) Fluent v6 é usado para modelar o

comportamento da superfície livre.

Mendes (2007) modelou numericamente o escoamento com superfície livre

ocasionado pelo movimento harmônico lateral ou vertical de um reservatório parcialmente

cheio de água. Duas possíveis condições iniciais para a superfície livre são definidas: repouso

absoluto ou uma onda senoidal estacionária. A agitação não-linear da superfície livre é

simulada primeiramente através do método da integral de contorno desenvolvido por Dold &

Peregrine (1986) e adaptado para as condições de contorno descritas acima. Os resultados

numéricos são validados por meio de soluções analíticas clássicas da literatura, tais como

Moiseyev (1958) e Chester (1968). Características interessantes da superfície livre são obtidas

e discutidas, tais como a formação de ondas de Faraday e o efeito bazuca (Longuet-Higgins

2001b).

Em uma segunda etapa, com a motivação de se explorar os seus recursos e

compará-lo com os estudos feitos na primeira parte, o código comercial ANSYS CFX é

empregado ao problema de sloshing horizontal e vertical, supondo um modelo homogêneo

para a superfície livre. Neste caso, as simulações são comparadas com os experimentos

realizados por Bredmose et al. (2003). Uma excelente concordância entre os resultados

numéricos e experimentais é encontrada.

Nicolato (2009) modelou numericamente o escoamento com superfície livre

ocasionado pelo movimento harmônico e não-harmônico lateral horizontal de um reservatório

parcialmente cheio de água – 1480 x 480 x 750 mm3 (comprimento x largura x altura). A

altura inicial da interface da água foi definida como 155 mm. O tanque é inicialmente

representando por uma caixa bidimensional impermeável, com um fundo plano horizontal e

duas paredes verticais. A condição de simetria é aplicada para as paredes normais à direção do

32

movimento do fluido. A condição inicial para a superfície livre de repouso absoluto é

definida.

O código comercial de fluidodinâmica computacional CFD, ANSYS CFX é

empregado ao problema de sloshing, aplicando-se um modelo multifásico homogêneo. A

malha computacional espacial é refinada próxima a superfície livre e próxima às paredes dos

reservatórios com o objetivo de melhor prever o perfil não-linear do movimento do fluido e

facilitar a convergência.

Os resultados numéricos são validados por meio de soluções analíticas clássicas

da literatura, tais como Moiseyev (1958) e Chester (1968). Em uma segunda etapa, o código

comercial ANSYS CFX é empregado ao problema de sloshing horizontal. Neste caso, as

simulações são comparadas com os experimentos realizados por Bredmose et al. (2003). Uma

excelente concordância entre os resultados numéricos e experimentais é encontrada e um

estudo dos parâmetros de simulação é realizado. O modelo é então extrapolado para demais

aplicações, como a modificação da geometria do reservatório conforme recomendações usuais

apresentadas em handbooks de construção offshore, de forma a minimizar o impacto das

ondas nas paredes dos reservatórios.

2.4 Método de volumes finitos

O Método de Volumes Finitos (MVF) consiste na discretização do domínio fluido

em volumes de controle. Determina-se uma distribuição de pontos, denominados pontos

nodais, dentro da geometria de estudo de maneira que cada nó é envolto por um volume de

controle, ou célula. É comum determinar volumes de controle próximos às extremidades do

domínio, de forma que os limites físicos coincidam com os limites dos volumes de controle.

A solução do problema de escoamento (velocidade, pressão, fração volumétrica

do fluido) é definida no centróide de cada célula. A precisão da solução é governada pelo

número de células do domínio, além da marcha temporal definida para os problemas

transientes. Com o aumento do número de nós, aumenta-se também a demanda por recursos

computacionais para a solução do problema. Para obtenção de uma boa precisão, portanto, é

necessário conseqüentemente, um grande tempo de processamento.

Denomina-se como malha, ou grade, a distribuição destas células no domínio.

Uma malha otimizada, geralmente, não é uniforme, caracterizada por maiores refinamentos

onde há grandes variações das propriedades entre os pontos nodais. Quanto à organização

33

desses elementos, classifica-se a malha como estruturada e não-estruturada. Quanto à forma

dos volumes de controle, é usual a divisão em elementos tetraédricos e hexaédricos.

Pode-se associar resumidamente o método dos volumes finitos aos seguintes

passos:

• Integração das equações de conservação em todas as células do domínio de

interesse.

• Substituição de diferentes esquemas de aproximação por diferenças finitas em

termos da equação integrada representando os processos do escoamento, como

convecção, difusão e fontes. Este processo converte as equações integrais em um

sistema de equações algébricas.

• Solução das equações algébricas por um método iterativo.

O primeiro passo distingue o método dos volumes finitos das demais técnicas de

CFD, o que acarreta na conservação das principais propriedades para cada elemento,

relevando uma das principais vantagens do esquema. A conservação de uma propriedade

genérica de um escoamento, ψ, por exemplo, um componente do campo de velocidades,

dentro do volume de controle pode ser expresso como o balanço entre os vários fatores que

tendem a contribuir ou reduzir a propriedade.

O segundo passo se refere à determinação dos esquemas adequados de

aproximação das derivadas presentes nas equações aplicáveis aos volumes de controle por

termos lineares. A seleção destes modelos de aproximação linear deve ser feita levando-se em

consideração as características do fenômeno de transporte envolvido. A seguir são listadas

algumas aproximações utilizadas:

• Aproximação por diferenças centrais, recomendado para problemas de caráter

difusivo, em que o valor da propriedade ψ é tomado levando-se em consideração a

contribuição das células adjacentes, sem considerar o fluxo do escoamento;

• Esquemas upwind e exponencial, aplicáveis principalmente a problemas

advectivos, em que o valor de uma propriedade genérica ψ é calculado

considerando-se a contribuição do nó vizinho de maior influência analisando-se a

direção do fluxo;

34

• Esquema upstream, que é uma formulação intermediária, podendo se adequar ao

caráter difusivo ou advectivo, aproximando-se do esquema de diferenças centrais

ou upwind.

Quanto à resolução das equações integradas no domínio temporal, podem-se

definir os esquemas de interpolação como as formulações totalmente implícita, implícita e

explícita (vide figura 2.4), representados pela equação geral a seguir:

0)1( Pψϕϕψψ ϕ −+= (2.7)

  W  P  E 

W  P  E t 

t+Δt W P E

W P Et

t+ΔtW P  E 

W P  E t 

t+Δt

  Explícita (φ=0)  Totalmente Implícita (φ=1)  Implícita (0<φ<1)

Figura 2.4 – Esquemas de interpolação temporal de uma propriedade genérica ψ

Normalmente o campo de velocidades é desconhecido, e os seus valores são

determinados pela solução global do problema junto com todas as demais variáveis do

escoamento. A solução das equações de transporte nas direções x, y e z e da equação da

continuidade apresenta alguns obstáculos, entre eles:

• Os termos advectivos da equação de momentum contêm termos não-lineares;

• Todas as equações são intrinsecamente acopladas, já que cada componente da

velocidade aparece nas equações de conservação de momentum e da continuidade.

Caso o escoamento considerado seja compressível, a equação da continuidade

pode ser usada como uma equação de transporte para a massa específica e temperatura através

do uso de uma equação de estado p=p(ρ,T). Quando o escoamento considerado é

incompressível, a massa específica é constante e não depende da pressão. Neste caso o

35

acoplamento entre a velocidade e a pressão introduz uma restrição na solução do campo de

escoamento: se o campo de pressões correto é aplicado às equações da continuidade, o campo

de velocidades resultante deve satisfazer a continuidade.

O acoplamento entre o campo de velocidades e a pressão, junto com as

não-linearidades presentes nas equações de momentum, pode ser resolvido iterativamente

através de alguns algoritmos como:

• SIMPLE (Semi-Implicit Method for Pressure Linked Equation) – Patankar &

Spalding (1972);

• SIMPLER (SIMPLE Revised) – Patankar (1980);

• SIMPLEC (SIMPLE Consistent) – Vand Doormal & Raithbay (1984);

• PISO (Pressure Implicit with Splitting of Operators) – ISSA (1986).

Para a solução desses algoritmos, costuma-se utilizar um recurso de discretização

conhecido como arranjo colocalizado para os componentes da velocidade e pressão. A ideia

consiste em avaliar as variáveis escalares, como a pressão, massa específica, temperatura, nos

pontos nodais, mas calcular os componentes da velocidade nas faces de cada célula.

2.5 ANSYS CFX

O programa de fluidodinâmica computacional ANSYS CFX apresenta as

seguintes etapas para resolver o problema de valor inicial:

• Pré-processamento;

• Solução;

• Pós-processamento.

Na fase de pré-processamento é realizada a modelagem do problema, inserindo

dados referentes à geometria do reservatório, características dos materiais, propriedades

termodinâmicas, características de transporte, de processo, do fluido e do sólido. Ainda

durante esta fase é realizado a discretização da malha, a atribuição de parâmetros numéricos

referentes a velocidade e precisão dos resultados, os modelos de cálculo para turbulência,

36

temperatura e pressão. A segunda fase (solução) é onde todas as características descritas no

pré-processamento se interagem, executando-se o solver do software. No pós-processamento,

a visualização dos resultados é realizada.

Com o programa ANSYS CFX é possível realizar a modelagem em um ambiente

multifásico, ou seja, pode-se realizar simulações de um ambiente onde se encontra presente

mais de um fluido. Para o problema proposto por este trabalho os dois fluidos simulados no

reservatório são a água e o ar. O ANSYS CFX representa as diferentes fases do ambiente

pelas letras gregas minúsculas α, β, γ etc.

Em geral, uma quantidade subscrita com α, β, γ, etc, se refere a um valor de uma

fase particular. Por exemplo, a fração volumétrica de α é denotado com αφ . Assim, o volume

αV ocupado pela fase α é um pequeno volume V ao redor de um ponto da fração de volume

αφ dado por:

VV αα φ= (2.8)

O número total de fases da simulação é definido como Np. A fração de volume de

cada fase é denotada de αφ .

O efeito da massa específica é então definido como:

ρφρ αα = (2.9)

Esta é a massa atual por unidade de volume de fase α, dado que a fase α somente

ocupa uma fração de volume, ou seja, a massa de α por unidade de grandeza do volume do

fluido.

A massa específica da mistura é dada por:

∑= ααφρρm (2.10)

A pressão total em uma simulação multifásica é definida como:

37

∑+=α

ααα ρφ 2

21 UPP stattot (2.11)

Esta definição é usada para ambos os escoamentos compressível e incompressível.

Num escoamento multifásico Euleriano dois diferentes sub-modelos são viáveis: o

modelo homogêneo e o modelo de transferência inter-fluido (não homogêneo).

No modelo não homogêneo o momento de transferência interfacial, calor e massa

são diretamente dependentes da superfície de contato entre as duas fases. Isto é caracterizado

pela área da interface por unidade de volume entre a fase a e fase b, sabendo-se que a

densidade de área superficial é Aαβ.

A transferência interfacial pode ser modelada usando ora partículas ou modelos

misturados. Esse promove diferentes descrições algébricas para a densidade da área

interfacial.

O modelo de partícula da transferência interfacial entre duas fases assume que

uma das fases é contínua (fase α) e outra é dispersa (fase β). A área de superfície por unidade

de volume é então calculada considerando que a fase β esteja presente como partícula esférica

de diâmetro dβ. Usando este modelo, a área de contato da interface é:

β

βαβ

φd

A6

= (2.12)

Os coeficientes de transferência adimensionais da interfase podem ser associados

aos termos da partícula do número de Reynolds e do número de Prandtl. Estes são definidos

usando o diâmetro da partícula e as propriedades da fase contínua, como a seguir:

α

βαβααβ μ

ρ dUU −=Re (2.13)

38

α

αααβ

μkCp

=Pr (2.14)

onde μα, Cpα e kα são as viscosidade, capacidade calorífica específica e a condutividade

térmica da fase contínua α.

No modelo homogêneo, um campo comum de escoamento é distribuído em todos

fluidos, como também outros campos relevantes tais como: temperatura e velocidade. Isto

permite fazer algumas simplificações no modelo multiflúido resultante do modelo

homogêneo.

Para um dado processo de transporte, o modelo homogêneo assume que a

quantidade transportada (com exceção da fração volumétrica) seja semelhante em todas as

fases, isto é:

ϕφα = 1≤ α ≤ NP (2.15)

Devido à quantidade transportada estar separada no escoamento homogêneo

multifásico, é suficiente resolver os campos comuns usando equações de transporte de volume

do que resolver individualmente equações de transporte da fase.

As equações de transporte volumétricas podem ser deduzidas somando as

equações individuais de transporte de fase, dando uma única equação do transporte para φ:

( ) ( ) SUt

=∇Γ−∇+∂∂ φφρρφ . (2.16)

onde

∑=

=Np

1ααα ρφρ , ∑

=

=Np

UU1

1α

ααα ρφρ

, ∑=

Γ=ΓNp

1αααφ (2.17)

39

O modelo homogêneo não necessita ser aplicado em todas as equações, por

exemplo, a velocidade do campo pode ser modelada como não homogênea, mas acoplado

com um modelo de turbulência homogênea. Por outro lado, um campo de velocidade

homogênea pode ser acoplado com um campo de temperatura não homogêneo. Variáveis

adicionais homogêneas são também permitidas no ANSYS CFX.

3 Modelagem e Controle

3.1 Introdução

A teoria de controle é a ferramenta da engenharia que, aplicada a um sistema

dinâmico (planta), possível de ser modelado matematicamente, permite a implementação de

sub-sistemas que levem uma ou várias variáveis do sistema a um estado predeterminado em

um tempo finito (Ogata, 2003). Duas estratégias podem ser identificadas: controle em malha

aberta e controle em malha fechada. Por exemplo, tratando-se de um sistema com uma

entrada e uma saída (single input single output – SISO), o controlador em malha aberta não

faz uso de qualquer informação sobre a saída. Já no caso de um controlador em malha

fechada, um sensor providencia informação sobre a saída (realimentação) para ser comparado

com um valor de referência, de tal modo que o erro entre estes sinais é processado pelo

controlador que finalmente fornecerá um sinal de atuação na entrada do sistema a ser

controlado. Na Fig. 3.1 se ilustram estes conceitos.

Figura 3.1 – Diagrama de blocos de um sistema de controle (a) em malha aberta, (b) em malha fechada.

Com o avanço da eletrônica e principalmente da aparição dos microprocessadores,

microcomputadores e microcontroladores, foi necessário levar em conta um novo parâmetro,

a saber, a freqüência de amostragem. Trata-se, neste caso, de lidar com sinais que devem

Saída Referência

- +

Ganho Planta Saída

Referência

(a)

(b)

Realimentação

u(t) e(t) Controlador Planta

40

passar do domínio contínuo para o domínio discreto ou digital, pois tem que se apresentar no

formato entendível pelos microprocessadores. Nesse domínio é necessário mudar a ferramenta

matemática da transformada de Laplace para a transformada Z, sendo que os procedimentos

de análise e projeto de sistemas de controle são, em geral, equivalentes. Entretanto, algumas

características são únicas nos sistemas de controle digital. A Fig. 3.2 ilustra um diagrama de

blocos de um sistema típico de controle digital em malha fechada, onde se nota a presença de

dois blocos característicos: o conversor analógico-digital (A/D) e o conversor digital-

analógico (D/A), os quais servem de interface entre os domínios dos sinais envolvidos no

sistema de controle.

Figura 3.2 – Diagrama de blocos de um sistema de controle digital em malha fechada.

3.2 Controladores

Diversos tipos de controladores podem ser encontrados na literatura. Aqueles

projetados inicialmente no domínio contínuo podem ser discretizados seguindo alguns dos

procedimentos já consagrados na literatura. Contudo, é possível projetar controladores

estritamente no domínio discreto, com características únicas, como por exemplo o controlador

digital dead-beat. Por outro lado, dependendo da planta, a estratégia de controle pode se

referir a um sistema SISO (single input – single output) ou a um sistema MIMO (multiple

input - multiple output). Neste último caso, uma matriz de realimentação de estados é

projetada ao invés de uma função de transferência para o controlador. Como a planta

abordada em este trabalho pode ser tratada de forma simplificada como se fosse um sistema

de primeira ordem em série com um integrador, então este item se limitará a descrever

brevemente os controladores básicos para sistemas SISO, os quais continuam sendo

encontrados em muitas aplicações industriais. Desta forma, podemos mencionar os seguintes:

• Controlador do tipo Proporcional (P);

Saída Referência

- +

Controlador Planta

Realimentação A/D

D/A

MICROCONTROLADOR

41

• Controlador do tipo Integral (I);

• Controlador do tipo Proporcional-Integral (PI);

• Controlador do tipo Proporcional-Derivativo (PD);

• Controlador do tipo Proporcional-Integral-Derivativo (PID).

3.2.1 Controlador proporcional (P)

É aquele onde o sinal da sua saída é diretamente proporcional ao sinal de erro

entre a referência e o sinal realimentado. Desta forma, considerando a Fig. 3.1 (b), este pode

ser representado matematicamente como:

.e(t)ku(t) p= (3.1)

onde kp é o ganho do controlador proporcional

3.2.2 Controlador integral (I)

É aquele onde a derivada temporal do sinal de saída é diretamente proporcional ao

sinal de erro conectado a sua entrada. Dessa forma este controlador pode ser

matematicamente representado por:

.e(t)kdt

du(t)i= (3.2)

Integrando a equação anterior, esta pode ser reescrita como:

∫=t

0i e(t)dt.ku(t) (3.3)

onde ki é o ganho do controlador integral

42

3.2.3 Controlador proporcional-integral (PI)

É aquele onde a saída do controlador resulta da combinação dos controladores

proporcional e integral descritos acima. Dessa forma esse controlador pode ser

matematicamente representado por:

∫+=t

oi

pp e(t)dt

Tk

.e(t)ku(t) . (3.4)

onde ki = kp / Ti e Ti é denominado de tempo integral.

3.2.4 Controlador proporcional-derivativo (PD)

É aquele que resulta da combinação de um controlador proporcional e outro cuja

saída é diretamente proporcional à derivada do erro. A representação matemática fica dada

por:

dtde(t).Tk.e(t)ku(t) dpp += . (3.5)

onde kd = kp . Td e Td é denominado de tempo derivativo.

3.2.5 Controlador proporcional-integral-derivativo (PID)

É aquele controlador composto pela combinação dos controladores proporcional,

integral e derivativo. É um dos controladores mais utilizados na indústria, freqüentemente

utilizado quando se deseja uma resposta transitória rápida e um erro de regime nulo. Desta

forma, pode ser representado como:

dtde(t).Tke(t)dt

Tk

.e(t)ku(t) dp

t

oi

pp ++= ∫ (3.6)

43

A versão discreta de cada um destes controladores pode ser obtida com auxilio de

alguns algoritmos entre os quais se podem mencionar: feedback e Tustin (Barcsack, 1995).

3.3 Acionamento com Estratégia de Controle PWM

A saída dos controladores acima mencionados são sinais ideais. Numa

implementação prática, este sinal não é conectado diretamente à planta, uma vez que, em

geral, os níveis de potência dos circuitos que implementam os controladores são diferentes

daqueles que acionam efetivamente a planta. Por exemplo, tratando-se de um motor, como é o

caso neste trabalho, a sua fonte de tensão deve proporcionar valores de corrente da ordem de

Ampères enquanto os circuitos eletrônicos somente lidam com correntes da ordem de décimos

de Ampère. Deste modo, os sinais de saída dos controladores somente podem servir de

referência para o acionamento do motor.

Em relação aos acionamentos, estes podem estar baseados em circuitos lineares,

utilizando transistores operando na região linear, ou em circuitos chaveados, onde os

transistores utilizados operam nas regiões de corte e saturação. A vantagem destes últimos

reside no fato de que o consumo de potência nas chaves se reduz drasticamente se comparado

com os acionamentos lineares. Por esta razão, neste trabalho foi utilizado um circuito de

potência chaveado.

Na figura 3.3 se mostra esquematicamente o arranjo do circuito de potência que

aciona o motor utilizado neste trabalho. Nele se podem apreciar quatro chaves semicondutoras

do tipo IGBT(isolated gate bipolar transistor), que são comandadas para se comportar como

interruptores fechados ou abertos. Este arranjo é chamado de ponte completa ou simplesmente

ponte H. A estratégia que comanda as chaves, ou estratégia de chaveamento, é feita de tal

forma que nunca se permita que duas chaves no mesmo braço fiquem em curto-circuito

(IGBT1 e IGBT2 ou IGBT3 e IGBT4). Quando as chaves IGBT1 e IGBT4 estão fechadas, as

outras devem estar abertas e vice-versa.

A seqüência de aberturas e fechamentos pode ser implementada de diversas

formas (Rashid, 1993) sendo que a utilizada neste trabalho é denominada de modulação por

largura de pulsos ou PWM, do inglês: Pulse Width Modulation. Nessa estratégia tem-se como

referência um período fixo de tempo T para o qual, a soma dos instantes de tempo em que

uma chave fica fechada (Ton) e aberta (Toff) é igual a esse período.

Desta forma, o motor alterna entre valores positivos e negativos de tensão, como

ilustrado na Figura 3.4, onde T1 é o instante em que as chaves IGBT1 e IGBT4 permanecem

fechadas enquanto IGBT2 e IGBT3 permanecem abertas. Durante o tempo T2 os estados das

44

chaves se invertem. Neste caso a tensão de alimentação Vcc é igual a 50V e o período de

referência vale 4 segundos ou, em outras palavras, a frequencia de alternância vale 0,25Hz. O

sentido de rotação e a velocidade do motor dependem do valor médio da tensão imposta ao

motor.

Figura 3.3 – Acionamento do motor em ponte completa.

Tensão no motor x Tempo

-60

-40

-20

0

20

40

60

0 1 2 3 4 5 6 7 Tempo (s)

Tensão (V)

Figura 3.4 – Tensão no motor x Tempo com uma estratégia PWM.

motor

IGBT3 IGBT4

IGBT1 IGBT2

Vcc Terra

T1 T2

T

45

A relação percentual entre o tempo em que uma chave semicondutora fica fechada

(T1) e o período T é chamada de ciclo de serviço D. Desta forma, a tensão média no motor

depende do ciclo de serviço. Por exemplo, quando D=50% a tensão média é zero e o motor

fica sem movimento. Se referindo às chaves IGBT1 e IGBT4, para ciclos de serviço maiores

do que 50% o sentido de rotação do motor é oposto ao sentido que acontece com um ciclo de

serviço menor do que 50%. A resolução do ciclo de serviço é definida pelo microcontrolador

ou microprocessador sendo utilizado.

3.4 Modelo da planta

Desde o ponto de vista da teoria de controle, o objeto a ser controlado é

denominado de planta e o seu modelo matemático deve oferecer suficiente informação para

permitir o projeto do controlador. Neste trabalho, o modelo da planta se refere a uma

plataforma de ferro acionada, através de uma polia e um parafuso sem fim, por um motor de

corrente contínua governado por um circuito de potência em ponte H como descrito no item

3.3. Em função das não linearidades da estrutura mecânica devido aos atritos presentes no

parafuso sem fim e no próprio eixo do motor, foi adotado um modelo experimental seguindo

o procedimento de Queiroz em 2007.

Para isto, o motor foi submetido a um degrau de tensão e foi observada a resposta

da sua velocidade, medida com auxilio de um codificador incremental. Uma resposta,

normalizada, similar à mostrada na Figura 3.5 foi observada. Esta é semelhante à que teria um

sistema de primeira ordem, o que significa que para atingir 95% do seu valor final o sistema

demora aproximadamente 3 vezes o valor da sua constante de tempo Tp.

Desta forma a função de transferência da planta G(s) que relaciona a velocidade

do motor com a tensão aplicada ao mesmo pode ser aproximada como:

1.sTkG(s)

p += , (3.7)

onde k é a relação, em regime permanente, entre a velocidade da mesa e a tensão aplicada e s

é a variável de Laplace.

46

Figura 3.5 – Resposta ao degrau unitário de um sistema de primeira ordem.

Para determinar os parâmetros do modelo foi aplicada a máxima tensão ao motor,

conseguida com um ciclo de serviço igual a 100%, que corresponde a uma tensão de 86V.

Verificou-se que a mesa demorou 0,32 segundos para chegar a 95% da sua velocidade

máxima de 0,1675 m/s. Com isto 3Tp=0,32, de onde Tp = 0,106 segundos e k = 0.001861

(m/s)/V.

Logo:

10,106s0,001861G(s)

+= (3.8)

Ou também:

9,5s0,01848G(s)+

= (3.9)

Como a variável de interesse é a posição da mesa, então a função de transferência

que relaciona a tensão do motor com a posição da mesa, P(s), pode ser obtida multiplicando

G(s) pela função de transferência de um integrador (1/s). Chamando de U(s) à função de

transferência da tensão aplicada ao motor e de X(s) à função de transferência da posição da

mesa, então P(s) pode ser escrita como:

9,5)s(s0,01848

U(s)X(s)P(s)

+== (3.10)

0 0,2 0,4 0,6 0,8

1 1,2

0 1 2 3 4 5 6

tempo

0,95

3Tp

U(s) P(S)

X(s)

47

3.5 Projeto do controlador

Para o projeto de um controlador, C(s) na figura 3.6, algumas técnicas clássicas

podem ser encontradas na literatura, sendo que deve ser levada em conta a forma da

implementação. Neste trabalho, trata-se de uma implementação digital, pois foi utilizado um

computador e o ambiente de programação Labview (Labview 7, 2003).

Sendo assim, existem duas possibilidades. A primeira consiste em projetar o

controlador no domínio contínuo (variável de Laplace) e depois convertê-lo para um

controlador discreto. A segunda consiste em discretizar o modelo da planta e projetar o

controlador no domínio discreto (variável Z). Neste trabalho foi adotada a primeira estratégia.

Como o modelo aproximado do sistema é simples, a técnica de projeto adotada

foi baseada no denominado lugar das raízes (Philips, 1996). Embora este modelo inclua um

integrador e seja possível projetar um simples controlador proporcional, de tal modo que após

ser aplicada uma referência em degrau o erro de regime seja nulo, outros controladores foram

avaliados em vista de que as referências de posição devem ter uma característica alternada

para gerar ondas no reservatório apoiado na mesa móvel. Além disso, outra razão para esta

avaliação reside na não linearidade por efeito do atrito mecânico que não foi contemplado no

modelo.

A seguir, será feita uma apresentação dos possíveis controladores a ser projetados,

onde, para facilitar este procedimento, a função de transferência do sensor de posição, H(s),

foi considerada como sendo um ganho unitário.

Figura 3.6 – Diagrama de blocos de um sistema de controle em malha fechada.

3.5.1 Controlador proporcional

Na figura 3.7 se observa o lugar das raízes do sistema de controle de posição com

um controlador proporcional. Este se origina nos pólos da planta e tende para valores

complexos conjugados com parte imaginária que tende para infinito. Neste caso, com

X(s) Referência

- +

H(s)

E(s) U(S) C(s) P(s)

48

qualquer ganho proporcional kp, Eq. 3.1, o sistema será estável e a velocidade de resposta do

controlador dependerá da posição dos pólos escolhidos.

-20 -15 -10 -5 0-1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

Real Axis

Imag

Axi

s

Figura 3.7 Lugar das raízes do sistema de controle com ganho proporcional.

3.5.2 Controlador proporcional derivativo

Neste caso, aplicando a transformada de Laplace à Eq. 3.5, a função de

transferência do controlador assume a forma: )1/T(sTKC(s) ddp += . Desta forma o lugar das

raízes incluirá um zero, cuja posição dependerá da escolha de Td. Na Figura 3.8 mostram-se

os lugares das raízes, com características estáveis, para os casos em que este zero fique à

direita do pólo -9,5 (a) e à esquerda do mesmo (b). No primeiro caso, a resposta do sistema

estará limitada pela presença de um pólo (em malha fechada) na região que fica entre o zero e

o pólo na origem. Já no segundo caso, a rapidez do sistema poderá aumentar

significativamente uma vez que o lugar das raízes se desloca à esquerda do pólo -9,5.

-20 -15 -10 -5 0-1

-0.5

0

0.5

1

Real Axis

Imag

Axi

s

-20 -15 -10 -5 0

-6

-4

-2

0

2

4

6

Real Axis

Imag

Axi

s

(a) (b)

Figura 3.8 Lugar das raízes do sistema com controlador proporcional derivativo.

49

Uma variante desse controlador consiste em adicionar um filtro passa baixa, cujo

pólo seja maior que o pólo da planta, com o objetivo de atenuar a amplificação de eventuais

ruídos de alta freqüência que acompanhem o sensor de posição. Na figura 3.9 se mostra a

mudança do lugar das raízes quando se inclui um filtro com pólo em -30rad/s e outro em -

40rad/s.

-50 -40 -30 -20 -10 0-30

-20

-10

0

10

20

30

Real Axis

Imag

Axi

s

-50 -40 -30 -20 -10 0-30

-20

-10

0

10

20

30

Real Axis

Imag

Axi

s

(a) (b)

Figura 3.9 Lugar das raízes do sistema com controlador proporcional derivativo e filtro.

3.5.3 Controlador proporcional integral

Aplicando a transformada de Laplace à Eq. 3.4, a função de transferência do

controlador assume a forma: )/s1/TKp(sC(s) i+= . Desta forma, o lugar das raízes incluirá

um pólo na origem e um zero, cuja posição dependerá da escolha de Ti. Na Figura 3.10

mostram-se os lugares das raízes, com características estáveis, para os casos em que este zero

fique à direita do pólo -9,5 (a) e (b). Nestes casos o sistema será estável e a resposta do

sistema será comparável à obtida com um controlador proporcional (Fig. 3.7). No caso (c), o

zero fica à esquerda do pólo, sendo que o sistema será instável pois o lugar das raízes se

desloca à direita do eixo imaginário.

50

-15 -10 -5 0-15

-10

-5

0

5

10

15

Real Axis

Imag

Axi

s

-15 -10 -5 0-10

-5

0

5

10

Real Axis

Imag

Axi

s

(a) (b)

-15 -10 -5 0-15

-10

-5

0

5

10

15

Real Axis

Imag

Axi

s

(c)

Figura 3.10 Lugar das raízes do sistema com controlador proporcional integral.

3.5.4 Controlador proporcional integral derivativo

Aplicando a transformada de Laplace à Eq. 3.6, a função de transferência do

controlador assume a forma: b)/ssaK(sC(s) 2 ++= , onde K, a e b dependem dos

parâmetros do controlador. Desta forma, o lugar das raízes incluirá um pólo na origem e dois

zeros, cujas posições no lugar das raízes dependeram da escolha dos parâmetros do

controlador. Supondo que estes ficassem no eixo real, três possibilidades existem. Os lugares

das raízes são mostrados na figura 3.11. Na figura 3.11 (a) os zeros ficam à direita do pólo da

planta -9,5. Na figura (b) um zero fica à direita e o outro à esquerda desse pólo e na figura (c)

ambos zeros ficam à esquerda daquele pólo. No primeiro caso, com um ganho apropriado, a

velocidade da resposta será comparável à obtida com o controlador proporcional. No segundo

caso será comparável à obtida com o controlador proporcional derivativo e, no terceiro caso,

51

apesar do lugar das raízes se deslocar à esquerda do pólo da planta,, haverá um pólo

dominante entre este pólo e o zero mais próximo dele, limitando sua velocidade à obtida com

o controlador proporcional derivativo.

-10 -5 0-2

-1.5

-1

-0.5

0

0.5

1

1.5

2

Real Axis

Imag

Axi

s

-30 -20 -10 0

-15

-10

-5

0

5

10

15

Real Axis

Imag

Axi

s

(a) (b)

-80 -60 -40 -20 0-40

-30

-20

-10

0

10

20

30

40

Real Axis

Imag

Axi

s

(c) Figura 3.11 Lugar das raízes do sistema com controlador proporcional integral derivativo.

A partir das observações dos possíveis candidatos a controlador pode-se concluir

que estes poderão ser o controlador proporcional, o controlador proporcional integral e o

controlador proporcional derivativo.

52

3.5.5 Discretização do controlador contínuo

Como foi mencionado no inicio, a implementação do controlador da posição da

bancada de testes foi feita via computador. Neste caso, a aquisição da posição e o algoritmo

de controle são efetuados de forma síncrona, tendo como referência uma freqüência fs fixa ou,

em outras palavras, um período Ts fixo.

Entre as diversas técnicas de discretização do controlador contínuo, uma consiste

em modelar o controlador em série com um grampeador de ordem zero e em seguida, através

de tabelas encontrar a correspondente transformada Z, para finalmente encontrar a

correspondente equação de diferenças implementada no computador. Entretanto algumas

técnicas aproximadas também podem ser utilizadas com sucesso, como por exemplo, as

aproximações: forward (Euler), backward e trapezoidal, também chamada de aproximação de

Tustin (Astrom, 1984). Estas aproximações consistem na substituição da variável s, na função

de transferência do controlador, por uma função da variável z, descritas a seguir:

Método de Euler: sT1-zs =

Método Backward: sT21-zs =

Método Tustin: 1z1-z

T2s

s +=

Neste trabalho, a aproximação de Tustin foi adotada uma vez que esta tem a

vantagem de transformar o semiplano esquerdo de Laplace no circulo unitário do plano z e

facilita a análise de estabilidade.

Desta forma, com k sendo o número da iteração utilizada na implementação real,

as versões discretas dos controladores candidatos e as suas respectivas equações de

diferencias, terão as seguintes formas: Controlador Versão Contínua Versão Discreta Equação de diferência Proporcional KpC(s) = KpC(z) = u(k) = Kp * e(k)

Proporcional Integral

s)1/T(sKpC(s) i+

= 1)-(zB)(zAC(z) +

= u(k)=u(k-1) + A*e(k) +

AB*e(k-1)

Proporcional Derivativo )1/T(sTKC(s) ddp += 1)(zD)(zCC(z)

++

= u(k)= - u(k-1) + C*e(k) + CD*e(k-1)

Onde A = Kp (1 + Ts/(2* Ti)); B = (Ts – 2Ti)/(Ts + 2Ti); C = Kp (1 + 2 Td/Ts); D = (Ts–2Td)/(Ts+2Ts/Td).

53

Controle digital com LABVIEW

Reservatório

Mesa móvel

Eletrônica de Potência

Câmera

Motor CC Interface USB de Aquisição de dados e

Controle

x

4 Descrição da bancada de testes

A Figura 4.1 mostra uma ilustração simplificada da bancada de testes

experimentais montada, neste trabalho, para a geração de ondas em um reservatório. A seguir

serão descritos os seus componentes e também o programa de controle desenvolvido para os

experimentos.

Figura 4.1. Bancada de testes.

4.1 Componentes

Para este trabalho foi utilizada uma mesa móvel disponível num laboratório de

automação do departamento de engenharia mecânica da UFF. Esta consiste de uma estrutura

em aço com dimensões aproximadas de 1,5m x 1,0 m com possibilidades de controle em dois

graus de liberdade através de dois motores de corrente contínua. Entretanto, apenas um grau

de liberdade foi utilizado para a geração de ondas, como mostra a Fig. 4.2.

O acoplamento da mesa com cada motor é feito através de polias, um parafuso

sem fim e um mecanismo de deslocamento por esferas (Queiroz, 2007). Os motores trabalham

com uma tensão máxima de 90Vcc e possuem uma potência nominal de 1,2 HP.

54

Figura 4.2 – Mesa móvel e reservatório.

Dois reservatórios retangulares, de vidro, de dimensões diferentes foram

utilizados nos experimentos. Estes ficaram apoiados na mesa tendo-se o cuidado de selar as

partes superiores dos mesmos, para proteger os demais componentes. Os experimentos se

diferenciam pelo tipo de controle: em malha aberta ou em malha fechada, e também pelo tipo

de controlador.

Para acionar o motor acoplado à mesa, foi construída uma placa eletrônica de

potência que consiste de semi-drivers do fabricante Semikron e transistores do tipo IGBT

montados na configuração em ponte completa. Estes são comandados por sinais binários

gerados a partir de uma placa eletrônica auxiliar que converte um sinal analógico, fornecido

pelo conversor D/A do computador, para um sinal do tipo PWM (Pulse Width Modulation).

Uma placa eletrônica do fabricante National Instruments, modelo USB-6001, foi

utilizada como interface entre o computador e o acionamento eletrônico do motor. O software

Labview foi utilizado para a implementação dos algoritmos de controle de posição.

Para medir a posição da mesa, e permitir o controle em malha fechada, foi

utilizado um sensor do tipo “encoder incremental”. Maiores detalhes sobre este dispositivo

podem ser encontrados em (Queiroz, 2007).

Um sensor da marca Signal Work, modelo 112A22, foi fixado a 10 mm da parte

inferior da parede esquerda do reservatório para medir as pressões de impacto das ondas.

Conectado a este sensor, foi utilizado um módulo de amplificação do sinal do sensor, um

computador e o programa de aquisição e processamento de dados CATMAN.

Para a filmagem e fotografias do movimento das ondas geradas pelo experimento

foi utilizada uma câmera Nikon D90 e uma filmadora digital Sony com capacidade de 33

quadros por segundo.

55

4.2 Programas de controle da mesa

Para este trabalho, dois programas de controle foram desenvolvidos no ambiente

de desenvolvimento Labview. Este se caracteriza por apresentar duas janelas de trabalho. Na

primeira, que será denominada de interface gráfica, são visualizados botões e ícones que

permitem ao usuário mudar comandos de operação, tais como ligar ou desligar, ou ainda

alterar um sinal de referência. A segunda janela está associada ao algoritmo de controle. Nela

reside efetivamente o programa de controle, implementado através da conexão de diversos

ícones que implementam operações lógicas e aritméticas. Nesta janela também se encontram

blocos associados à inicialização dos periféricos do computador e à freqüência de

amostragem. A esta janela se lhe denominou de janela de ambiente de desenvolvimento.

O primeiro programa implementado se refere a um sistema de controle de

movimento em malha aberta. Neste caso, o programa envia um sinal de referência alternado,

para o conversor D/A da placa USB-6001, que é traduzido pela placa auxiliar em um ciclo de

serviço que alterna, de forma senoidal, em torno de 50%. Com isto a tensão elétrica no motor

alterna entre valores positivos e negativos que resultam em um movimento alternado da mesa

em torno de uma posição de equilíbrio. Estes experimentos foram denominados de N1 e N4.

Nas figuras 4.3 e 4.4 são apresentadas as janelas de ambiente de desenvolvimento

e de interface gráfica respectivamente, para este experimento. Na primeira se aprecia um

comando, “controle de tensão manual”, que está associado, na segunda janela, a um

potenciômetro de cabeça redonda com escala de 0 a 3. Também, à esquerda da primeira

janela, pode ser identificado o gerador de sinal alternado, cuja amplitude é somada a um valor

constante de 1,97V devido às características da placa auxiliar. Observa-se à direita da

primeira figura o bloco “DAQ Assistant”, onde são definidas as configurações da placa de

aquisição de dados. Os outros ícones servem para habilitar a visualização gráfica, da

referência de tensão senoidal do motor e a posição da mesa, e também da visualização

numérica de outros parâmetros.

56

Figura 4.3 – Ambiente de desenvolvimento

Figura 4.4 - Interface gráfica do programa.

O segundo programa desenvolvido se refere a um sistema de controle de

movimento em malha fechada, denominado de experimento N3. Este foi concebido para lidar

com os diversos tipos de controladores.

Na figura 4.5 é apresentada a janela de interface gráfica e nas figuras 4.6 e 4.7 são

apresentadas janelas ligadas ao ambiente de desenvolvimento, para este experimento.

57

Figura 4.5 - Interface gráfica do programa.

Figura 4.6 – Ambiente de desenvolvimento

58

Figura 4.7 – Ambiente de desenvolvimento (cont.)

Para este programa apresentamos aqui os principais módulos desenvolvidos.

O primeiro módulo é o responsável pela geração sinal senoidal de referência,

amortecido durante o primeiro tempo t1 ajustado pelo programa. Esse módulo basicamente

consiste de uma rampa temporizada multiplicado por um sinal de seno também ajustado para

terminar no instante de tempo t2 ajustado pelo programa.

Um segundo módulo importante desse programa é o módulo de controle, onde o

sinal de referência é comparado com a posição real da mesa, obtendo o erro entre os dois

sinais. Esse erro é enviado ao controlador do programa, que conforme discutido no capítulo 3

pode ser de vários tipos.

Um terceiro bloco está constituído por vários módulos que são os responsáveis

pelo envio e recebimento de sinais para a placa USB da National Instruments. Outros

módulos, não menos importantes, são: o módulo de contagem, desenvolvido em Fortran para

funcionar como o contador de pulsos que identifica em que posição se encontra a mesa; o

módulo de saturação, também desenvolvido em Fortran, responsável por limitar a saída de

tensão do controlador para níveis aceitáveis de disparo do ciclo de serviço do PWM e da placa

de interface; e por último a interface gráfica do programa, responsável pelo comando,

monitoramento e controle da mesa.

59

4.3 Experimento N3: teste em malha fechada

No intuito de se testar o controle em malha fechada, estabeleceu-se o seguinte

experimento – denominado por N3 – na bancada de testes. Foi utilizada como referência um

sinal de excitação do tipo senoidal com amplitude atenuada nos primeiros 6,0s (vide figura

4.8). Após os 6,0s, o movimento segue a função,

x(t) = 15 sen (3,142 t). (5.2)

Neste caso, a amplitude escolhida é de 15mm e a freqüência de 3,142rad/s ou

0,50Hz, com um período de 2,0s.

Posição de referência

-20,00

-15,00

-10,00

-5,00

0,00

5,00

10,00

15,00

20,00

02,

44,

87,

199,

59 1214

,416

,819

,221

,6 2426

,428

,831

,233

,535

,938

,340

,743

,145

,547

,950

,352

,755

,157

,559

,962

,3

Tempo (s)

Posiç

ão (m

m)

Figura 4.8. Referência de deslocamento da mesa para o experimento N3.

Para este experimento, o reservatório utilizado possui as seguintes dimensões

internas: 800 x 100 x 400 mm³ (comprimento, largura e altura). Vários testes foram realizados

nesta categoria, os quais foram denominados H1, H2, H3, H4, H5 e H6. Entretanto os

melhores resultados foram os testes H4 e H5. A diferença principal entre estes é que em H4 o

reservatório foi preenchido com 60mm de água, enquanto que em H5 o reservatório foi

preenchido com 70mm de água.

Foram desenvolvidos para estes experimentos, controladores do tipo proporcional

e proporcional-integral, entretanto ficou experimentalmente comprovado que o controlador

proporcional gerava menos vibrações na superfície livre do fluido, coerentes com as vibrações

observadas na própria mesa, como será mostrado mais adiante. Como a experiência com o

controlador PD foi semelhante à obtida com o controlador PI, os resultados com o controlador

60

PD não serão apresentados nesse trabalho. Deve-se lembrar que no projeto teórico dos

controladores, apresentado no capítulo 3, tinha-se assumido que o ganho do sensor era

unitário mas na implementação real o ganho equivalente do sensor foi de 1000. Com isto, os

ganhos utilizados no ambiente Labview se reduziram na proporção.

Inicialmente será apresentado o comportamento da posição da mesa utilizando

controladores PI. Três casos são apresentados, diferenciados pela amplitude dos ganhos. Em

ambos os casos os ganhos proporcional e integral foram iguais, o que significou a introdução

de um zero em –1. Os pólos de malha fechada dependem dos ganhos escolhidos. No caso A,

Kp = Ki = 1,6. No caso B, Kp = Ki = 0,79 e, no caso C, Kp = Ki = 0,3.

Na figura 4.9 é mostrado o lugar das raízes resultante e a posição dos pólos para

os casos A e B. A resposta dos controladores é mostrada nas figuras 4.10, 4.11 e 4.12, onde a

referência de posição cresce linearmente durante 4 períodos para logo ficar em um patamar

constante. Como percebido, o melhor comportamento da posição foi no caso C mas sempre

acompanhado de vibração embora menor que nos outros casos. Outros casos também foram

testados, verificando que a posição do zero do controlador tinha que ser à direita do pólo –9,5

da planta.

Depois da experiência com controladores PI, foram feitos testes experimentais

para diversos controladores proporcionais. Verificou-se que os ganhos em torno de 0,1 para

referências de posição alternadas com freqüências de 1 e 1,1 Hz foram apropriados, como

mostrado nas figuras 4.13 e 4.14.

(a) (b)

Figura 4.9. Lugar das raízes e pólos em malha fechada.

61

Figura 4.10. Resposta experimental para o caso A, onde a linha branca corresponde a

saída e a linha vermelha corresponde a referência de posição da mesa.

Figura 4.11. Resposta experimental para o caso B, onde a linha branca corresponde a

saída e a linha vermelha corresponde a referência de posição da mesa.

62

Figura 4.12. Resposta experimental para o caso C, onde a linha branca corresponde a

saída e a linha vermelha corresponde a referência de posição da mesa.

Figura 4.13. Resposta do controlador proporcional com ganho Kp = 0,1 e referência de

posição com freqüência de 1,1Hz.

63

Figura 4.14. Resposta do controlador proporcional com ganho Kp = 0,1 e referência de

posição com freqüência de 1,0Hz.

Como observado, com controladores proporcionais as vibrações na posição se

reduziram significativamente mas em compensação o erro de regime entre a referência e o

deslocamento real da mesa não foi nulo. Apesar desta característica se considerou preferível

movimentar o reservatório com o controlador proporcional.

Após a realização de vários testes experimentais e simulações optou-se por adotar

como ganho do controlador proporcional o valor de 0,07. Os resultados do comportamento

cinemático das ondas formadas, para os casos H4 e H5, são apresentados nas figuras 4.15 e

4.16, respectivamente.

64

Figura 4.15. Três instantes em seqüência do movimento horizontal H4 do reservatório.

t=28,7s

t=29,9s

t=29,2s

65

Figura 4.16. Três instantes em seqüência do movimento horizontal H5 do reservatório.

t=24,4s

t=25,7s

t=27,5s

66

Em relação à pressão medida na parede esquerda do reservatório, a figura 4.17

mostra o resultado para o experimento H4. Conforme pode ser observado, existe uma

seqüência de impactos na parede do reservatório compatível com o movimento da mesa e com

o período do sinal de excitação utilizado. Para os instantes de tempo entre 22 e 60 segundos

parece existir um certo padrão na medição da pressão. Entretanto, as variações do perfil

podem ser explicadas pela resolução do sensor e as vibrações na mesa. A figura 4.18 mostra o

resultado para o caso H5, sendo que parece existir um padrão entre os instantes 30 e 60

segundos.

-6,00E-01

-4,00E-01

-2,00E-01

0,00E+00

2,00E-01

4,00E-01

6,00E-01

8,00E-01

0 20 40 60 80 100 120

Tempo (s)

Pres

são

(psi

)

Figura 4.17. Pressão x tempo na parede do reservatório para o teste H4.

-3,00E-01

-2,00E-01

-1,00E-01

0,00E+00

1,00E-01

2,00E-01

3,00E-01

4,00E-01

5,00E-01

0 20 40 60 80 100 120

Tempo (s)

Pres

são

(psi

)

Figura 4.18. Pressão x tempo na parede do reservatório para o teste H5.

67

Com o objetivo de se avaliar as limitações do sensor de pressão e/ou do sistema

de aquisição de dados e de se verificar a capacidade de registro de um duplo pico de pressão,

característica do fenômeno de sloshing, foi testado manualmente um movimento lateral na

direção x, de forma semelhante ao experimento H10 de Bredmose et al. (2003), de tal forma

que permitisse a geração de ondas energéticas que impactariam em uma das paredes do

reservatório. A figura 4.19 ilustra o perfil de pressão encontrado, comprovando a capacidade

de detecção do duplo pico de pressão quando do run up na parede vertical.

Convém registrar que as pressões registradas durante os primeiros 20s nas figuras

4.17 e 4.18, e durante os primeiros 8s na figura 4.19, devem ser descartados por conta de um

atraso inerente do sensor de pressão utilizado.

-2,00E-01

-1,00E-01

0,00E+00

1,00E-01

2,00E-01

3,00E-01

4,00E-01

5,00E-01

0 2 4 6 8 10 12 14 16

Tempo (s)

Pres

são

(psi

)

Figura 4.19. Pressão x tempo na parede do reservatório para o teste manual.

68

5 Resultados

Os resultados experimentais e numéricos encontrados são apresentados e

discutidos nas seções a seguir. A condição inicial de repouso absoluto é empregada na

superfície livre, com o reservatório sendo agitado lateralmente com determinadas frequências

e amplitudes de oscilação. Em todos os casos, apenas movimentos horizontais na direção x

são considerados em nosso estudo.

Com o objetivo de se formar um benchmark sobre o assunto, primeiramente são

reproduzidas as simulações numéricas realizadas por Mendes (2007) e Nicolato (2009), que

compararam os experimentos H10, realizados por Bredmose et al. (2003), e resultados

numéricos obtidos por intermédio do código comercial ANSYS CFX. Utilizando-se o aparato

experimental relatado nos capítulos 3 e 4, dois estudos de caso – batizados de N1 e N4 – são

apresentados e discutidos. Os experimentos N1 e N4 são então comparados às

correspondentes soluções numéricas encontradas via ANSYS CFX.

5.1 Validação via experimento H10 (Bredmose et al. 2003)

Nesta seção, os resultados experimentais e numéricos de Bredmose et al. (2003)

para o fenômeno de sloshing horizontal (experimento H10) são comparados com soluções

numéricas obtidas através do código comercial ANSYS CFX, considerando-se um modelo

multifásico homogêneo para a descrição do movimento da superfície livre. Para maiores

detalhes do método numérico empregado e do software ANSYS CFX, veja o capítulo 2. A

figura 5.1 mostra a malha hexaédrica de 9207 elementos hexaédricos (18800 nós) utilizada

nessa simulação.

Nas paredes verticais e no fundo do reservatório foi aplicada a condição de não-

deslizamento, enquanto na região superior a condição de opening – livre transferência de

massa – foi empregada. O regime adotado foi o transiente, com passo de tempo de 0,001s,

para um tempo total de 9,0s, dentro da faixa do experimento de Bredmose et al. O

reservatório é modelado com dimensões 1,48 x 0,40 x 0,75 m3 e é parcialmente cheio de água,

com uma profundidade média de 0,155 m. A figura 5.2 ilustra o sinal utilizado nas simulações

numéricas. Os resultados para a superfície livre são mostrados em um referencial inercial.

69

Figura 5.1. Discretização empregada nas simulações com o software ANSYS CFX.

Figura 5.2. Movimento harmônico lateral da mesa empregado no experimento H10 de

Bredmose et al. (2003). O deslocamento horizontal dispx é dado em metros e o tempo t em

segundos.

A coluna à direita da figura 5.3 mostra uma sequência de 12 fotos extraídas do

experimento H10 (ver figura 11 de Bredmose et al.), enquanto as linhas tracejadas

representam as soluções numéricas utilizando-se um modelo de Boussinesq. À esquerda da

mesma figura, os resultados obtidos utilizando-se o código ANSYS CFX são apresentados.

As imagens cobrem o movimento da superfície livre resultante do “empurrão” ocasionado

pela parede à direita do reservatório, decorridos 6,21s do início do experimento (vide figura

5.2).

A onda gerada é captada até atingir a parede à esquerda, seguida do run up e

posterior reflexão. De uma forma geral, a concordância entre os resultados experimentais e

numéricos é excelente. A amplitude da onda gerada na parede vertical à direita decai até se

aproximar da parede oposta, quando se torna mais íngreme e inicia a sua ascensão. Note que

70

os resultados numéricos obtidos no CFX reproduzem mais fielmente os experimentos se

comparados aos perfis obtidos a partir do modelo de Boussinesq, notoriamente no momento

do run up.

Como reportado por Bredmose et al., o modelo de Boussinesq é falho em regiões

em que a batimetria não é suave, como, por exemplo, próximo à parede vertical. Uma outra

característica importante evidenciada nas simulações é a indicação, após o run up, de um jato

vertical em direção ao fundo do reservatório (ver figura 5.3, t = 8.04s), fato este não

reproduzido pelo modelo de Boussinesq.

Convém ressaltar ainda que os resultados numéricos apresentados não

contemplam os efeitos da tensão superficial, que segundo Jervis & Peregrine (1996) podem

acarretar em uma redução de até 10% da amplitude máxima em fenômenos não-lineares.

Nicolato (2009) considerou os efeitos da tensão superficial quando da simulação do

experimento H10 de Bredmose et al. Os perfis de onda numéricos encontrados (com e sem

tensão superficial) são praticamente idênticos, sendo que o caso com tensão superficial

curiosamente apresentou um pico ligeiramente maior, sem no entanto influenciar o perfil de

pressões na parede. Tal diferença pode ser atribuída à necessidade de uma melhor

discretização na região do run up, tendo em vista as pequenas escalas de comprimento

associadas aos fenômenos capilares.

Os resultados obtidos para o campo de pressões com o software ANSYS CFX são

comparados com as medições e cálculos efetuados por Bredmose et al., com uma

concordância excelente (ver figura 5.4). Percebe-se claramente a formação do tradicional

duplo pico de pressão quando da ocorrência do run up na parede vertical. Conforme já

discutido no capítulo 3 e 4, por limitações da bancada de testes, o experimento H10 não pode

ser reproduzido. Assim, dois estudos de caso foram estabelecidos visando a validação das

simulações numéricas, batizados por N1 e N4, que se seguem.

Maiores detalhes de como acontece o fenômeno do duplo pico de pressão podem

ser encontrados no capítulo 2, item 2.3 desse trabalho.

71

Figura 5.3. Comparação entre os resultados numéricos (obtidos a partir do código comercial

ANSYS CFX) e experimentais (experimento H10 de Bredmose et al. 2003).

72

Figura 5.3 (cont.). Comparação entre os resultados numéricos (obtidos a partir do código

comercial ANSYS CFX) e experimentais (experimento H10 de Bredmose et al. 2003).

73

(a)

(b)

Figura 5.4. (a) Comparação entre as pressões medidas (- - - - -) e calculadas (_______)

por Bredmose et al. (2003) para o experimento H10. (b) Distribuição de pressões obtida

através do código comercial ANSYS CFX.

74

5.2 Experimento N1

Para o experimento N1, foi aplicado um movimento horizontal harmônico no

reservatório,

x(t) = 11 sen (5,236 t), (5.1)

onde a amplitude escolhida é de 11mm, a freqüência é de 5,236rad/s ou 0,83Hz, com um

período de 1,2s. A amplitude da tensão aplicada ao motor foi proporcional a esta referência

(vide figura 5.5).

-15

-10

-5

0

5

10

15

0,03 0,12 0,21 0,3 0,39 0,48 0,57 0,66 0,75 0,84 0,93 1,02 1,11 1,2 1,29 1,38 1,47

Tempo (s)

Amplitude (mm)

Figura 5.5. Referência de deslocamento da mesa para o experimento N1.

As dimensões internas do reservatório foram de 393 x 194 x 265mm³

(comprimento, largura e altura). Não foi realizada a medição de pressão nas paredes do

reservatório. O reservatório foi preenchido com uma coluna de 50mm de água. As imagens

capturadas a partir da filmadora digital, para três instantes de tempo diferentes, são mostradas

na figura 5.6.

Essas imagens correspondem a uma onda deslocando-se da esquerda para direita.

Após o impacto da onda na parede do reservatório ocorre um jato vertical caracterizando o

fenômeno run up. Em seguida, um movimento de queda livre é iniciado, formando uma

depressão na superfície livre, que posteriormente se desloca da direita para a esquerda.

75

Figura 5.6. Três instantes em seqüência do movimento horizontal do reservatório.

5.3 Experimento N4

Considerando as limitações da mesa (referidas à máxima tensão aplicada ao

motor) e a taxa de amostragem máxima da placa de aquisição de dados da National

Instruments, esse experimento foi realizado na tentativa de se obter ondas mais energéticas,

que evidenciassem o duplo pico de pressão presente no fenômeno de sloshing. O experimento

N4 se caracterizou por ser em malha aberta, utilizando o reservatório do experimento N3, com

uma coluna de água de 50mm.

Para o experimento N4, foi aplicado um movimento horizontal harmônico no

reservatório,

x(t) = 90 sen (1,885 t), (5.1)

t=2,85s

t=3,03s

t=2,55s

76

onde a amplitude escolhida é de 90mm, a freqüência é de 1,885rad/s ou 0,3Hz, com um

período de 3,3s. O sinal de excitação do sistema é mostrado na figura 5.7. A figura 5.8 mostra

os perfis das ondas formadas.

Figura 5.7. Interface gráfica do experimento N4.

77

Figura 5.8. Três instantes em seqüência do movimento do reservatório com o

experimento N4.

t=9,4s

t=9,9s

t=13,1s

78

5.4 Simulação numérica e validação via experimento N1

A figura 5.9 exibe a malha utilizada na simulação numérica do experimento N1.

Foi utilizado um maior refinamento na região próxima à superfície livre, com uma malha

hexaédrica com 18480 nós e 9047 elementos. O movimento horizontal foi imposto em todo o

domínio fluido. A simulação foi realizada em regime transiente, com passo de tempo de

0,005s para um tempo total de 12s.

Figura 5.9. Malha hexaédrica utilizada no experimento N1.

A figura 5.10 compara cinco instantes de tempo extraídos do experimento (coluna

à direita) com os resultados numéricos correspondentes obtidos utilizando-se o software

ANSYS CFX (coluna à esquerda). Uma boa concordância entre os resultados é encontrada.

Pode-se notar a formação de uma onda de baixa amplitude decorridos 0,51 s, que se choca

contra a parede esquerda do reservatório decorridos 1,20s sem causar, no entanto, um jato

vertical. No tempo 2,10s a onda já se encontra com maior amplitude, agregando

características não-lineares ao escoamento, atingindo a formação de run up em t=2,80s. Com

a queda da coluna d’água (t=3,00s), há a formação de uma depressão próxima à parede

vertical, evidenciada tanto no experimento quanto na simulação numérica.

79

Figura 5.10. À esquerda, resultados numéricos obtidos com o ANSYS CFX. À direita,

sequência de fotos extraída do experimento N1.

5.5 Simulação numérica e validação via experimento N4

Para esta nova geometria foi utilizada uma malha hexaédrica com cerca de 20000

nós e 10000 elementos, bem como um refinamento na região próxima à superfície livre. A

simulação foi realizada em regime transiente, com passo de tempo de 0,005s para um tempo

total de 60s. A comparação entre os perfis da superfície livre estão evidenciadas na figura

5.11. Observa-se a existência de uma boa concordância entre a cinemática das ondas

experimentais e numéricas. No instante t=7,6s, após o impacto na onda gerada na parede

direita do reservatório, é formada uma seqüência de pequenas ondas refletidas, registradas no

instante t=8,1s.

t=0,51s

t=1,20s

t=0,51s

t=1,20s

t=2,10s t=2,10s

t=3,00s t=3,00s

t=2,80s t=2,80s

80

Figura 5.11. À esquerda, resultados numéricos obtidos com o ANSYS CFX. À direita,

sequência de fotos extraída do experimento N4.

t=7,6s t=7,6s

t=7,8s

t=8,0s

t=8,1s

t=7,8s

t=8,0s

t=8,1s

t=8,4s

t=8,7s

t=8,4s

t=8,7s

81

Figura 5.11. (cont.) À esquerda, resultados numéricos obtidos com o ANSYS CFX. À

direita, sequência de fotos extraída do experimento N4.

Nos experimentos, a aquisição da pressão apresentou uma série de limitações. Foi

utilizado um sensor de pressão que tem o funcionamento limitado a sinais de baixa

freqüência, como é o caso dos sinais gerados pelos experimentos N3 e N4. Dessa forma, os

resultados da aquisição de cerca de 20 segundos iniciais precisam ser desprezados, pois esse

período retrata o tempo que o sensor de pressão leva para estabilizar a medição, conforme

discutido na seção 4.3.

Apesar de ter sido verificado manualmente que o sensor registra o duplo pico de

pressão quando da presença do fenômeno de sloshing (vide figura 4.19 após 8s), a bancada de

testes não consegue reproduzir movimentos que gerem ondas energéticas. Assim, pode-se

observar na figura 5.12a que existe um histórico de variação da pressão na parede do

reservatório no experimento N4, mas que, no entanto, diverge dos resultados numéricos

encontrados com o programa ANSYS CFX (vide figura 5.12b). Maiores investigações

experimentais são necessárias para se identificar a faixa de oscilações que é possível de ser

imposta na mesa e que gere perfis de pressão com magnitudes que sejam compatíveis com o

sensor de pressão utilizado.

t=9,3s

t=11,3s

t=9,3s

t=11,3s

82

-6,00E-01

-4,00E-01

-2,00E-01

0,00E+00

2,00E-01

4,00E-01

6,00E-01

8,00E-01

0 10 20 30 40 50 60 70

Tempo (s)

Pres

são

(psi

)

(a)

(b)

Figura 5.12. Experimento N4: variação da pressão (a) medida e (b) calculada via CFX.

Tempo (s)

Pres

são

(Pa)

83

6 Conclusões

A aplicação de métodos numéricos a escoamentos com superfície livre tem se

mostrado extremamente eficiente na previsão de fenômenos lineares e, sobretudo, não-

lineares. Simulações numéricas, no entanto, necessitam ser validadas por soluções analíticas

ou experimentos em laboratório. Métodos analíticos em geral apresentam a desvantagem de

serem aplicáveis apenas a problemas de valor de contorno simples, mas em compensação

trazem soluções elegantes e compactas. Experimentos em laboratório têm a grande vantagem

de poderem tratar a configuração real de determinado escoamento. Entretanto, podem

alcançar elevado custo e, muitas vezes, não são viáveis de serem realizados por questões de

segurança ou devido à dificuldade de se reproduzir as condições reais.

Em trabalhos anteriores realizados no Laboratório de Dinâmica dos Fluidos

Computacional (LabCFD) da Universidade Federal Fluminense (UFF), resultados numéricos

foram validados através da comparação com resultados analíticos e numéricos (Mendes 2007,

Nicolato 2009), além da comparação com experimentos realizados no exterior (Bredmose et

al. 2003). Neste sentido, este trabalho focou no desenvolvimento de uma bancada

experimental para validação das simulações numéricas na UFF, considerando o problema de

geração e propagação de ondas em um reservatório.

O experimento consistiu na adaptação de uma mesa xy para ser controlada por

meio de um circuito eletrônico de potência que acionou um motor de corrente contínua. Este,

por sua vez recebe um sinal do tipo PWM proveniente de um circuito eletrônico auxiliar.

Todo o conjunto foi controlado através de um programa desenvolvido na plataforma Labview

que interage, através da placa USB-6001 da National Instruments, com os circuitos acima e

um sensor de posição. Esta configuração mostrou-se adequada apesar das limitações do

hardware.

Os programas desenvolvidos em computador possibilitaram o controle da mesa

tanto em malha aberta quanto em malha fechada. O programa desenvolvido para o controle

em malha aberta possibilitou que a mesa trabalhasse com diversas freqüências e diversas

amplitudes do sinal de excitação do motor. O programa desenvolvido para o controle da

posição da mesa em malha fechada foi restrito a menores freqüências e amplitudes de posição

devido à velocidade de aquisição de dados da placa National Instruments (interface analógica

digital). Esse programa possibilitou a utilização de controladores diversos como proporcional,

proporcional-integral e proporcional-derivativo. A partir dos resultados experimentais se

84

verificou que são necessários outros estudos para melhorar a modelagem da planta e testar

outros controladores com melhor desempenho do que os que foram observados.

As ondas geradas experimentalmente em diferentes reservatórios retangulares

foram simuladas numericamente e comparadas, apresentando uma boa conformidade com

relação aos perfis da superfície livre e à fase das ondas propagadas. Entretanto, a distribuição

de pressões registrada no sensor ao longo da parede do reservatório apresentou discrepâncias

consideráveis quando comparada aos resultados numéricos, excetuando o experimento H10 de

Bredmose et al. (2003). Tal diferença pode ser atribuída à faixa dos sinais gerados pelos

experimentos N3 e N4, todos de baixa freqüência. Um outro ponto negativo a ser evidenciado

é que não se conseguiu com o aparato experimental reproduzir o movimento do experimento

H10 de Bredmose et al. (2003), tampouco registrar o duplo pico de pressão característico do

fenômeno de sloshing.

Como trabalhos futuros são sugeridas as seguintes ações:

• Substituição da placa de aquisição e transmissão de dados por outra com maior

velocidade de operação.

• Utilização de um computador dedicado e com sistema operacional em tempo real,

com o intuito de operar com sinais de referência de posição de maior freqüência.

• Utilização de outros sensores de pressão de forma a aumentar a faixa de operação

permitindo a detecção de sinais de baixas freqüências.

• Implementação de outros graus de liberdade na mesa de trabalho.

• Substituição do sensor de posição do tipo encoder por outro de tipo absoluto.

• Aplicação de outros métodos numéricos ao problema de valor de contorno, tais

como elementos finitos e de contorno.

• Implementação do experimento para verificação de ondas ressonantes, conforme

trabalhos de Chestek (1968) e Moiseyev (1958). Foram calculadas as seguintes

freqüências de ressonâncias e amplitudes para número de onda igual a 1:

experimento N1 (f = 1,62 Hz, A = 1,23 cm), experimento N4 (f = 0,85 Hz, A =

2,34 cm).

85

7 Referências bibliográficas

ASTROM and WITTENMARK, Computer Controlled Systems, 1984.

CZESLAU L. BARCZAK – Controle digital de sistemas dinâmicos – Editora Edgard

Blucher, 1995

BREDMOSE, H., BROCCHINI, M., PEREGRINE, D.H. and THAIS, L., “Experimental

investigation and numerical modelling of steep forced water waves.” J. Fluid Mech., vol.

490, pp.217-249, 2003.

FALTINSEN, O.M. and TIMOKHA, A.N., “Sloshing.” Cambridge University Press,

2009.

IBRAHIM, R.A., “Liquid Sloshing Dynamics: Theory and Applications.” Cambridge

University Press. 948p, 2005.

KIM, Y., “Experimental and numerical analyses of sloshing flows.” J. Eng. Math., v.58,

pp.191-210, 2007.

LABVIEW 7 Express, National Instruments. USA, 2003.

LAMB, H., “Hydrodynamics.” Cambridge University Press, 1932.

MOREIRA, R.M., MENDES, A.A.O. and BACCHI, R.D.A., “The generation and

propagation of nonlinear waves in a reservoir.” Proceedings of the 11th Brazilian

Congress of Thermal Sciences and Engineering. Curitiba. pp.1-8, 2006.

NICOLATO, P.C. and MOREIRA, R.M., “Numerical Modeling of Water Wave Impact

on Reservoirs.” Proceedings of the 20th International Congress of Mechanical

Engineering. Gramado. pp.1-9, 2009.

OGATA, K., “Engenharia de Controle Moderno.” Editora Prentice Hall do Brasil, 2003.

86

PEREGRINE, D.H., “Water Wave Impact on Walls.” Ann. Rev. Fluid Mech., v.35, pp.1-

22, 2004.

PHILIPS, C.L., “Sistemas de Controle e Realimentação.” Makron Books, 1996.

QUEIROZ, M.T., “Automação de Mesa Coordenada Destinada a Máquina Fresadora

CNC.” Dissertação de Mestrado. UFF, 2007.

RASHID, M.H., “Power Eletronics Circuits, Devices, and Applications.” Prentice-Hall,

Inc., 1993.

TOCCI, R.J., “Sistemas Digitais.” Ed. Pearson Prentice Hall, 2003.

ABRAMSON, H. N., Slosh Suppression, NASA SP-8031, 1969.

ANDERSON, J. G.; SEMERCIGIL, S. E.; TURAN, O. F., A standing-wave-type

sloshing absorber to control transient oscillations. Journal of Sound and Vibration, v. 232,

n. 5, p. 839-856, 2000.

ANSYS CFX – Solver, Release 11.0: Theory; 2006 ANSYS Europe Ltd.

ATTARI, N.K.A.; ROFOOEI, F.R., On Lateral response of structures containing a

cylindrical liquid tank under the effect of fluid/structure resonances, Journal of Sound and

Vibration, v. 318, pp.1154-1179, 2003.

BATCHELOR, G. K., An Introduction to Fluid Dynamics. Cambrigde University Press,

2000.

CHEN, Y. G.; DJIDJELI, K.; PRICE, W. G., Numerical simulation of liquid sloshing

phenomena in partially filled containers. J. Computers & Fluids, v. 38, n. 1, p. 830-842.

Chester, W., 1968, Proc. Royal Soc, A306, 5, 2009.

CHO, J. R.; LEE, H. W.; HA, S. Y., Finite element analysis of resonant sloshing response

in 2-D baffled tank. Journal of Sound and Vibration, v. 288, p. 829-845, 2005.

87

CHUN, M. S.; KIM, M. H.; KIM, W. S.; KIM, S. H.; LEE, J. M., Experimental

investigation on the impact behavior of membrane-type LNG carrier insulation system.

Journal of Loss Prevention in the Process Industries, v. 22, p. 901-907, 2009.

COOKER, M,J.; PEREGRINE, D.H., Wave Impact Pressure and Effect Upon Bodies

Lying on the Bed, Coastal Engng, v.18, pp.205-229, 1992.

COX E. A. ; GLEESON J. P. ; MORTELL M. P., Nonlinear sloshing and passage

through resonance in a shallow water tank. Zeitschrift fur angewandte Mathematik und

Physik, v. 56, n. 4, p-645-480, 2005.

CRAIG, K. J.; KINGSLEY, T. C., Design optimization of containers for sloshing and

impact. Journal of Structural and Multidisciplinary Optimization, v. 33, n. 1, p. 71-87,

2007.

GAVRILYUK, I.; LUKOVSKY, I.; TROTSENKO, Y.; TIMOKHA, A. Sloshing in a

vertical circular cylindrical tank with an annular baffle. Part 1. Linear fundamental

solutions. Journal of Engineering Mathematics, v. 54, p. 71-88, 2006.

GODDERIDGE, B.; TURNOCK, S.; TAN, M.; EARLS, D. An investigation of

multiphase CFD modelling of lateral sloshing tank. Computers and Fluids, v. 38, n. 2, p.

183-193, 2009.

GRACZYK, M.; MOAN, T. A probabilistic assessment of design sloshing pressure time

histories in LNG tanks. Journal of Ocean Engineering, v. 35, p. -834-855, 2008.

JERVIS, M.; PEREGRINE, D. H. Overtopping of Waves at Sea Wall: A Theorical

Approach. In: PROC. 25TH. CONF. COASTAL ENGNG, 2., ASCE. p. 2192 - 2205.

1996.

KYOUNG, J. H.; HONG, S. Y.; KIM, J.W.; KWANG, J. B. Finite-element

computational of wave impact load due to a violent sloshing. J. of Ocean Engineering, v.

32, n.1, p. 2020-2039, 2005.

88

LEE, D. H.; KIM, M. H.; KWON, S. H.; KIM, J. W.; LEE, Y. B. A parametric sensitivity

study on LNG tank sloshing loads by numerical simulations. Journal of Ocean

Engineering, v. 34, n. 1, p. 3-9, 2007.

MA, Q. W.; YAN, S. Quasi ALE finite element method for nonlinear water waves.

Journal of Computational Physics, v. 212, n. 1, p. 52-72, 2006.

MARINTEK. MARINTEK and the safe and efficient LNG sea transport. Disponível em

<http://www.sintef.no/Home/Marine/MARINTEK/MARINTEK-

Publications/MARINTEK-Review-No-1---April---2004/MARINTEK-and-the-Safe-and-

Efficient-LNG-Sea-Transport/>. Acesso em 03 set. 2009.

MARINTEK. Sloshing and Impact Loads in Membrane Tanks. Disponível em

<http://www.sintef.no/Home/Marine/MARINTEK/MARINTEK-

Publications/MARINTEK-Review-No-1---April---2004/Sloshing-and-Impact-Loads-in-

Membrane-Tanks/>. Acesso em 03 set. 2009.

MENDES, André A. O. Geração de ondas não-lineares em um reservatório. Trabalho de

Projeto Final, 2º semestre de 2006, Universidade Federal Fluminense.

MOREIRA, R.M. & PEREGRINE, D.H. Nonlinear interactions between a free surface

flow with surface tension and a submerged cylinder. Journal of Fluid Mechanics (to

appear), 2009.

MOISEEV N. N. Some problems in the hydrodynamics of surface waves, in Advances in

Mechanics in USSR over 50 years, Moscow, Nauka, v. 2, p. 55–72, 1970.

MOISEEV N. N.; RUMNYANTSEV V. V. Dynamic Stability of Bodies Containing

Fluid, Applied Physics and Eng./An. Int. Series, 6, H. N. Abramson, ed., Springer-Verlag,

1968.

89

PETROBRAS. Tipos de Plataformas. Disponível em

<http://www2.petrobras.com.br/Petrobras/portugues/plataforma/pla_tipo_plataforma.htm

> Acesso em: 20 out. 2009.

MOISEYEV, N.N. On the theory of nonlinear vibrations of a liquid. Prikl. Mat. Mekh.

22, pp. 612–621, 1958.

RABIER, S.; MEDALE, M. Computation of free surface flows with a projection FEM in

a moving mesh framework. J. Comput. Methods Appl. Mech. Engng., v. 192, n. 1, p.

4703-4721, 2003.

ROYON-LEBBEAUD, A.; HOPFINGER, E. J.; CARTELLIER, A. Liquid sloshing and

wave breaking in circular and square-base cylindrical containers. J. Fluid Mechanics, v.

577, p. 467-494, 2007.

VERSTEEG, H.K.; MALALASEKERA, W. An introduction to computational fluid

dynamics - The finite volume method. EUA: Longman Scientific & Technical, 1995.

WEI, G.; KIRBY, J. T. A Time Dependent Code for Extended Boussinesq Equations. J.

Waterway Port Coastal Ocean Engng, v. 121, p. 251 – 261. 1995.

WEI, G.; KIRBY, J. T.; GRILLI, S.T.; SUBRUMANYA, R. A Fully Nonlinear

Boussinesq Model for Surface Waves. Journal Fluid Mechanics, v. 294, p. 71-92. 1995.

YONGHWANM, K. Experimental and numerical analyses of sloshing flows. Journal of

Engineering mathematics, v. 58, n. 1-4, p. 191-210, 2007.

YUELONG, Z.; XING, J. A boundary element method for small-amplitude viscous fluid

sloshing in couple with structural vibration. Acta mechanica solida sinica, v. 12, n. 1, p.

9-21, 1999.

ABRAMSON, H.N. The dynamic behaviour of liquids in moving containers. NASA

SP106. 1966.

90

ANDERSON, A.; DIVER, D. A., PEREGRINE, D. H. Slender, Steep and Breaking

Water Waves. In: PROC. 5TH INTL WORKSHOP ON WATER WAVES AND

FLOATING BODIES (ed. P. A. Martin), pp. 8-13. Manchester University, 1990.

BENJAMIN, T. B.; URSELL, F. The Stability of Plane Free Surface of a Liquid in

Vertical Periodic Motion. Proc. R. Soc. Lond, A 225, p.505-515. 1954.

BILLINGHAM, J. Nonlinear sloshing in zero gravity. Journal Fluid Mechanics, v.464,

p.365–391. 2002.

CHAPMAN, G.J.D.; PORTER, R. Second-Order Sloshing Over an Arbitrary Bed.

Journal Fluid Mechanics, v.524, p. 331- 355. 2005.

DAVIS, A. M. J. Two Dimensional Oscillations in a Canal of Arbitrary Cross Section.

Proc. Camb. Phil. Soc, v.61(1), p.827 – 846. 1965.

DAVIS, A. M. J. Short Surface Waves in a Canal; Dependence of Frequency on The

Curvatures and Their Derivatives. Q.J. Mech. Appl. Maths, v.27(1), p.523 – 535. 1974.

DOLD, J.W. An Efficient Surface-Integral Algorithm Applied to Unsteady Gravity

Waves. J. Comp. Phys., v.103, p.90-115. 1992.

DOLD, J.W.; PEREGRINE, D.H. An efficient boundary-integral method for steep

unsteady water waves. In Numer. Meth. Fluid Dyn. Ii, (Eds. K.W. Morton & M.J.

Baines), p. 671–679. 1986.

FALTINSEN, O.M. A Nonlinear Theory of Sloshing in Rectangular Tanks. J. Ship Res.,

v.18, p. 224 – 241. 1974.

GOODRIDGE, C.L.; HENTSCHEL, H.G. E.; LATHROP, D.P. Breaking Faraday

Waves: Critical Slowing of Droplet Ejection Rates. Phys. Rev. Lett., v. 82, p.3062 –

3065. 1999.

91

HATTORI, M.; ARAMI, A.; YUI, T. Wave Impact Pressure on Vertical Walls Under

Breaking Waves of Various Types. Coastal Engng, v.22, p.79-114. 1994.

JERVIS, M.; PEREGRINE, D. H. Overtopping of Waves at Sea Wall: A Theorical

Approach. In: PROC. 25TH. CONF. COASTAL ENGNG, 2., ASCE. p. 2192 - 2205.

1996.

JIANG, L.; PERLIN, M.; SCHULTZ, W.W. Period Tripling and Energy Dissipation of

Breaking Standing Waves. Journal Fluid Mechanics, v.369, p.273 - 299. 1998.

LONGUET-HIGGINS, M. S. Vertical Jets From Standing Waves. Proc. R. Soc. Lond., v.

A 457, p. 495-510. 2001a.

LONGUET-HIGGINS, M. S. Vertical Jets From Standing Waves; The Bazooka Effect.

IUTAM SYMPOSIUM ON FREE SURFACES FLOWS, p. 195 – 203. 2001b.

LONGUET-HIGGINS, M. S.; COCKELET, E.D. The Deformation of Steep Surface

Waves on Water. I. A Numerical Method of Computation. Phil. Trans. R. Soc. Lond., v.

A 350, p.1-26. 1976.

LONGUET-HIGGINS, M. S.; DOMMERMUTH, D.G. On The Breaking of Standing

Waves By Falling Jets. Phys. Fluids, v.13, p.1652 -1659. 2001a.

LONGUET-HIGGINS, M. S.; DOMMERMUTH, D.G. Vertical Jets From Standing

Waves II. Proc. R. Soc. Lond., v. A 457, p.2137-2149. 2001b.

LONGUET-HIGGINS, M. S.; DRAZEN, D.A. On Steep Gravity Waves Meeting a

Vertical Wall: A Triple Instability. Journal Fluid Mechanics, v.466, p.305 -318. 2002.

MCIVER, P.; PEREGRINE, D.H. Comparison of Numerical and Analytical Results for

Waves That Are Starting to Break. IN PROC. CONF. HYDRODYNAMICS IN OCEAN

ENGNG, TRONDHEIM, p. 203 – 215, 1981.

92

MERCER, G. N.; ROBERTS, A. J. Standing Waves In Deep Water: Their Stability and

Extreme Form. Phys. Fluids, v. A4, p.259-269. 1992.

MERCER, G. N.; ROBERTS, A. J. The Form Standing Waves on Finite Water Depth.

Wave Motion, v.19, p.233-244, 1994.

MILES, J.; HENDERSON, D. Parametrically Forced Surface Waves. Annu. Rev. Fluid

Mech., v. 22, p.143 – 165, 1990.

NIECKELE, A.A. Computação Avançada em Fenômenos de Transporte: Programa CFD

para Solução de Problemas Bidimensionais Elípticos. Rio de Janeiro, 1996.

OCKENDON, J.R.; OCKENDON, H. Resonant Surface Waves. Journal Fluid

Mechanics, v.59, p. 397- 413, 1973.

OCKENDON, H.; OCKENDON, J.R; JOHNSON, A.D. Resonant Sloshing in Shallow

Water. Journal Fluid Mechanics, v.167, p. 465- 479, 1986.

PACKHAM, B. A. Small-Amplitude Waves in a Straight Channel of Uniform Triangular

Cross-Section. Q.J. Mech. Appl. Maths , v.33, p.179 – 187, 1980.

STOKES, G.G. On Theory of Oscillatory Waves. Camb. Trans., V.8, p.441- 455. 1847.

SUCHITHRA, N.; KOOLA, P.M. A Study of Wave Impact on Horizontal Slabs. Ocean

Engng, v.22, p.687-697, 1995.

TELES DA SILVA, A.F., CHACALTANA, J.T.A. & MOREIRA, R.M. Nonlinear

sloshing of gravity waves in a reservoir. Proc. COBEM 2003, 2003, p.1–6.

TOPLISS, M. E. Water Wave Impact on Structures. Bristol, 1994. Tese (PhD

Dissertation) - School of Mathematics University of Bristol, UK, 1994.

URSELL, F. Proc. Roy. Soc., v. A 214, p.79, 1952.

93

WOOD, D.J.; PEREGRINE, D.H. Wave Impact Beneath a Horizontal Surface. 25th

INTERNATIONAL CONFERENCE COASTAL ENGINEERING. ASCE , p. 1-9, 1996.

WU, G.X.; EATOCK TAYLOR, R. Finite Element Analysis of Two-Dimensional

Nonlinear Transient Water Waves. Appl. Ocean Res., v.16, p.363 – 372, 1994.

WU, G.X.; EATOCK TAYLOR, R. Numerical Simulation of Sloshing Waves in a 3 D

Tank Based on a Finite Element Method. Appl. Ocean Res.,v.16, p.363 – 372, 1998.

ZHANG, S.; YUE, D.K.; TANIZAWA, K. Simulation of Plunging Wave Impact on a

Vertical Wall. J. Fluid Mech., v. 327, p.221 – 254, 1996.

ISSA, R.I..Solution of the Implicity Discretised Fluid Flow Equations by Operator-

Splitting,, J. Comput. Phys., Vol. 62, pp. 40-65, 1986.

PANTAKAR, S. V. Numerical Heat Transfer and Fluid Flow, Hemisphere Publising

Corporation, Taylor & Francis Group, New York, 1980.

PANTAKAR, S.V. and SPALDING, D.B. A Calculation Procedure for Heat, Mass and

Momentum Transfer in Three-dimensional Parabolic Flows, Int. J. Heat Mass Transfer,

Vol. 15, p. 1787, 1972.

VAN DOORMAL, J.P. and RAITHBY, G.D. Enhancements of the SIMPLE Method for

Predicting Incompressible Fluid Flows, Numer. Heat Transfer, Vol. 7, pp. 147-163, 1984.

94

8 Apêndice

Projeto do controlador PI

Abaixo apresentamos o desenvolvimento do controlador proporcional-integral

utilizado no experimento.

Dada a planta a ser controlada, definida pela expressão,

)5,9(01861,0

)()(

+=

sssHsX (A.1)

e dado um controlador PI definido pelas expressões abaixo:

szsk

skikps

skikp )1( +

⇒+

⇒+ (A.2)

Gostaríamos de calcular o lugar das raízes definindo, definindo inclusive os

ganhos kp e ki do controlador.

Considerando que a função de transferência em malha aberta do sistema é dada

pela expressão a seguir,

)5,9()1(01861,0.)( 2

'

++

=ss

zsksH (A.3)

Podemos dizer que esta expressão possui:

• Dois pólos em zero

• Um pólo em -9,5

• Um zero em Z1

95

Da função A.3 acima apresentada, podemos concluir que a expressão do

controlador, poderia ser definida como sendo:

ssksC )2.()( +

= (A.4)

Assim sendo, podemos proceder com o cálculo dos ganhos Kp e ki comparando as

duas expressões anteriores, da seguinte forma:

kiskpksk +=+ .2. (A.5)

kkpkki

== 2

ou kkp

kzki== .1

(A.6)

Utilizando a função RLOCFIND no programa Matlab, escolhemos como ganho

k1 para o controlador o valor de 39,8648.

Considerando ainda que conforme a equação A.5 acima,

01861,0.1 kk = (A.7)

Temos que:

12,214201861,08648,39

==k (A.8)

96

Ou seja, temos então que considerando Z1=2, as constantes proporcionais e

integrais do controlador são:

24,4284

12,2142

=

=

ki

kp (A.9)

Dessa forma, um possível controlador para o sistema pode ser definido como

sendo:

sssC )2.(12,2142)( +

= (A.10)

Discretização do controlador contínuo

“A discretização é o processo pelo se pode converter sinais contínuos em

discretos. Diz-se que um sistema controlado digitalmente é discreto ou discretizado. Quando o

sistema e as variáveis envolvidas são contínuas, deve-se convertê-las em discretas por meio

da amostragem no tempo.” (Czeslau 1995).

Sabendo que o nosso controlador e sistema é contínuo no tempo, e que o

controlador da planta seria implementado por um software através de um computador,

utilizamos os procedimentos analíticos da transformação bilinear ou Método de Tustin para a

discretização do compensandor ou controlador digital da planta.

Dado o controlador desenvolvido no item anterior desse capítulo definido por:

sssC )2.(12,2142)( +

= (A.11)

Onde ainda kp=2142,12 e ki=4284,24, podemos escrever o controlador da

seguinte forma:

97

skikpsC +=)( (A.12)

ssC 24,428412,2142)( += (A.13)

Segundo a teoria da transformação bilinear ou Método de Tustin, basta substituir

“s” pela expressão A.14 que poderemos produzir um controlador discretizado D(z), com base

na fórmula de integração trapezoidal.

Dado:

112

+−

=zz

Ts (A.14)

Substituindo “s” na equação A.13 do controlador, teremos:

)1(2)1.(.24,428412,2142

)()()(

−+

+==z

zTzEzUzD (A.15)

Considerando agora que o período que o irá trabalhar é de 0,004 segundos, o

equivalente a uma freqüência de250 Hz, podemos escrever que:

)1(2)1.(14,17)1.(24,4284

)()()(

−++−

==z

zzzEzUzD (A.16)

)1(57,8.57,812,2142.12,2142

)()()(

−++−

==z

zzzEzUzD (A.17)

98

)1(55,2133.69,2150

)()()(

−−

==zz

zEzUzD (A.18)

Multiplicando o numerador e denominador por 1−z , temos:

1

1

1.55,213369,2150

)()()( −

−

−−

==z

zzEzUzD (A.19)

11 ).(.55,2133)(.69,2150).()( −− −=− zzEzEzzUzU (A.20)

Aplicando a transformada inversa de Z, chegamos a equação abaixo que é a

expressão que representa o controlador da planta discretizado:

)1(.55,2133)(.69,2150)1()( −−=−− kekekuku (A.21)

)1(.55,2133)(.69,2150)1()( −−+−= kekekuku (A.22)