( PDF ) Projeto de turbinas hidráulicas axiais com parametrização da geometria, equação de equilíbrio radial e técnicas de otimização

Download PDF

ads:

UNIVERSIDADE FEDERAL DE ITAJUBÁ

INSTITUTO DE ENGENHARIA MECÂNICA

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

DISSERTAÇÃO DE MESTRADO

Projeto de Turbinas Hidráulicas Axiais

com Parametrização da Geometria, Equação de

Equilíbrio Radial e Técnicas de Otimização

Autor: Rodrigo Barbosa da Fonseca e Albuquerque

Orientador: Prof. Dr. Nelson Manzanares Filho

Co-Orientador: Prof. Dr. Waldir de Oliveira

Itajubá, Agosto de 2006

ads:

Livros Grátis

http://www.livrosgratis.com.br

Milhares de livros grátis para download.

UNIVERSIDADE FEDERAL DE ITAJUBÁ

INSTITUTO DE ENGENHARIA MECÂNICA

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

DISSERTAÇÃO DE MESTRADO

Projeto de Turbinas Hidráulicas Axiais

com Parametrização da Geometria, Equação de

Equilíbrio Radial e Técnicas de Otimização

Autor: Rodrigo Barbosa da Fonseca e Albuquerque

Orientador: Prof. Dr. Nelson Manzanares Filho

Co-Orientador: Prof. Dr. Waldir de Oliveira

Curso: Mestrado em Engenharia Mecânica

Área de Concentração: Dinâmica dos Fluidos e Máquinas de Fluxo

Dissertação submetida ao Programa de Pós-Graduação em Engenharia Mecânica como

parte dos requisitos para obtenção do Título de Mestre em Engenharia Mecânica.

Itajubá, Agosto de 2006

M.G. – Brasil

ads:

UNIVERSIDADE FEDERAL DE ITAJUBÁ

INSTITUTO DE ENGENHARIA MECÂNICA

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

DISSERTAÇÃO DE MESTRADO

Projeto de Turbinas Hidráulicas Axiais

com Parametrização da Geometria, Equação de

Equilíbrio Radial e Técnicas de Otimização

Autor: Rodrigo Barbosa da Fonseca e Albuquerque

Orientador: Prof. Dr. Nelson Manzanares Filho

Co-Orientador: Prof. Dr. Waldir de Oliveira

Composição da Banca Examinadora:

Prof. Dr. João Roberto Barbosa - ITA

Prof. Dr. Ariosto Bretanha Jorge – IEM/UNIFEI

Prof. Dr. Waldir de Oliveira - IEM/UNIFEI

Prof. Dr. Nelson Manzanares Filho, Presidente - IEM/UNIFEI

Dedicatória

Ao tio Deco.

Agradecimentos

Mais uma vez, o Prof. Nelson conduziu um trabalho interessante e criativo. Além

disso, seu “espírito de Telê Santana”, jogando com categoria e pra frente, tornou bastante

agradável e produtiva sua orientação.

A contribuição do Prof. Waldir para este trabalho vai além dos diversos artigos

cedidos e das discussões de alguns tópicos. Sendo um legítimo entusiasta da área de máquinas

de fluxo, é também um grande educador. A ele deve-se boa parte do destacado desempenho

dos alunos da UNIFEI na área de turbomáquinas nos últimos anos.

Minha família sempre fez sacrifícios visando os estudos meu e de meus irmãos. Difícil

deve ser para quem não conta com seu carinho. Vó Badi e vô Orlando toleraram minha

presença em sua casa por oito anos.

Sem os amigos e colegas da Universidade, o esforço não teria a menor graça. E sem o

futebol no Colégio das Irmãs, o estresse e a preocupação fariam mais estragos do que calvície

ou caspa.

Jana, eu te amo!

“Quanto mais você compreende as limitações de seu conhecimento

acerca de um assunto, mais apto você está para superá-las”

J. D. Denton

“Concept is most important, not the details”

Theodore von Karman

Resumo

ALBUQUERQUE, R. B. F. (2006), Projeto de Turbinas Hidráulicas Axiais com

Parametrização da Geometria, Equação de Equilíbrio Radial e Técnicas de

Otimização, Itajubá, 94p. Dissertação (Mestrado em Dinâmica dos Fluidos e Máquinas

de Fluxo) - Instituto de Engenharia Mecânica, Universidade Federal de Itajubá.

Este trabalho apresenta o desenvolvimento de uma metodologia computacional de baixo

custo para o projeto otimizado de turbinas hidráulicas axiais. A metodologia foi desenvolvida

usando-se um modelo de escoamento quase-bidimensional, com correlações empíricas para as

perdas e desvios nas grades. O estudo baseou-se nos princípios de conservação de massa,

energia e quantidade de movimento, considerando-se a equação de equilíbrio radial para se

obter um campo de escoamento mais realista. Tendo-se em vista um número reduzido de

variáveis de projeto, os ângulos de montagem, as razões corda-passo e os arqueamentos das

pás foram parametrizados em termos de seus valores nas estações do cubo, meio e ponta. O

algoritmo de projeto otimizado – solver e método de otimização – foi implementado num

programa escrito na linguagem MatLab

. Esse programa busca geometrias básicas que

maximizam o rendimento da turbina, dadas a vazão, rotação, restrições para a altura de queda

e faixas para as variáveis de projeto.

Duas técnicas de otimização foram aplicadas: um método de busca local, baseado em

gradiente, e um algoritmo populacional. O método de gradiente é uma Programação

Quadrática Seqüencial padrão, usando-se a função

fmincon do MatLab

, que busca mínimos

locais partindo-se de uma estimativa inicial. Para as buscas globais, apresentam-se os

Algoritmos de Busca Aleatória Controlada, algumas versões e suas vantagens no problema de

projeto otimizado.

Um exemplo de aplicação da metodologia é apresentado e discutido. Trata-se da

otimização de uma turbina hélice tubular existente, previamente ensaiada em bancada de

testes. Soluções ótimas são comparadas com o projeto original da turbina, mostrando-se as

melhorias de desempenho, segundo a modelagem hidrodinâmica. Sugestões para o

aperfeiçoamento da metodologia são dadas no final.

Palavras-chave

Máquina de Fluxo, Turbina Hidráulica Axial, Modelagem de Perdas e Desvio,

Parametrização da Geometria, Técnica de Otimização, Projeto Otimizado.

Abstract

ALBUQUERQUE, R. B. F. (2006), Design of Axial-Flow Hydraulic Turbines with

Geometry Parameterization, the Radial Equilibrium Equation and Optimization

Techniques, Itajubá, 94p. MSc. Dissertation - Instituto de Engenharia Mecânica,

Universidade Federal de Itajubá.

This work presents the development of a low cost computational methodology for the

conceptual design optimization of axial-flow hydraulic turbines (propeller turbines). The

methodology has been developed with a quasi-two dimensional flow model, employing

empirical correlations for cascade losses and flow deviations. The study is based on the

conservation principles for mass, energy and momentum. The radial equilibrium equation is

included in order to achieve a more realistic flow field. For reducing the number of design

variables, the runner blading stagger, chord-pitch ratio and camber are parameterized in terms

of their values at the hub, mean and tip stations. The design optimization algorithm has been

coded in MatLab™ language. This code searches for a basic geometry that maximizes the

turbine efficiency, given the design flow rate, rotational speed and bounds for the design

variables and also for the available head.

Two optimization techniques have been applied: a gradient based local search method

and a population set-based global search algorithm. The gradient based technique is a

standard Sequential Quadratic Programming, using the

fmincon function from MatLab™,

which searches for local minimizers starting from an initial point. For the global searches, it is

presented the Controlled Random Search Algorithms, some of their versions and their

advantages in the design optimization problem.

An application example of the methodology is presented and discussed for the

optimization of a real tube type propeller turbine, previously tested in a laboratory rig. The

optimized solutions are compared with the original turbine design, showing the performance

improvements, according to the hydrodynamic modeling. Suggestions for methodology

improvements are also made.

Keywords

Turbomachine, Axial-Flow Hydraulic Turbine, Loss and Deviation Modeling,

Geometry Parameterization, Optimization Technique, Design Optimization.

Sumário

SUMÁRIO_________________________________________________________________i

LISTA DE FIGURAS_______________________________________________________iii

LISTA DE TABELAS ______________________________________________________ v

SIMBOLOGIA ____________________________________________________________vi

LETRAS LATINAS ________________________________________________________vi

LETRAS GREGAS _______________________________________________________ vii

SUPERESCRITOS________________________________________________________viii

SUBSCRITOS____________________________________________________________viii

ABREVIATURAS _________________________________________________________ix

SIGLAS __________________________________________________________________ x

CAPÍTULO 1 _____________________________________________________________ 1

INTRODUÇÃO E MOTIVAÇÃO ____________________________________________ 1

CAPÍTULO 2 _____________________________________________________________ 5

FORMULAÇÃO DO PROBLEMA DE OTIMIZAÇÃO__________________________ 5

CAPÍTULO 3 _____________________________________________________________ 8

PARAMETRIZAÇÃO DA GEOMETRIA DAS PÁS ____________________________ 8

3.1 Modelo Geométrico ------------------------------------------------------------------------------8

3.2 Parametrização-----------------------------------------------------------------------------------15

CAPÍTULO 4 ____________________________________________________________ 20

CÁLCULO DO ESCOAMENTO EM UMA TURBINA HIDRÁULICA AXIAL_____ 20

4.1 Equação de Equilíbrio Radial------------------------------------------------------------------23

4.2 Equilíbrio Radial Após o Distribuidor --------------------------------------------------------27

4.3 Equilíbrio Radial Após o Rotor----------------------------------------------------------------31

4.4 Correlação de Desvio; Triângulos de Velocidade-------------------------------------------35

4.5 Correlações de Perdas; Rendimento da Turbina---------------------------------------------38

4.5.1 Perda no distribuidor ----------------------------------------------------------------------40

4.5.2 Perda no rotor ------------------------------------------------------------------------------40

4.5.3 Perda no tubo de sucção-------------------------------------------------------------------41

4.5.4 Perda mecânica-----------------------------------------------------------------------------42

4.5.5 Correlação de Soderberg------------------------------------------------------------------42

4.5.6 Integração das perdas e da potência absorvida pelas pás -----------------------------44

4.5.7 Rendimento e altura de energia da turbina ---------------------------------------------46

CAPÍTULO 5 ____________________________________________________________ 49

MÉTODOS DE OTIMIZAÇÃO_____________________________________________ 49

5.1 Introdução ----------------------------------------------------------------------------------------49

5.2 Programação Quadrática Seqüencial----------------------------------------------------------52

5.3 Algoritmos de Busca Aleatória Controlada (CRSA)----------------------------------------54

5.3.1 Algoritmo básico---------------------------------------------------------------------------55

5.3.2 Algumas versões ---------------------------------------------------------------------------55

CAPÍTULO 6 ____________________________________________________________ 65

RESULTADOS E DISCUSSÃO _____________________________________________ 65

6.1 Análise do Projeto Inicial de Souza (1989) --------------------------------------------------65

6.2 Resultados do SQP-Fmincon-------------------------------------------------------------------72

6.3 Resultados do CRSA-CRSI --------------------------------------------------------------------79

CAPÍTULO 7 ____________________________________________________________ 85

CONCLUSÕES E SUGESTÕES ____________________________________________ 85

7.1 Conclusões----------------------------------------------------------------------------------------85

7.2 Sugestões -----------------------------------------------------------------------------------------87

REFERÊNCIAS BIBLIOGRÁFICAS ________________________________________ 90

iii

Lista de Figuras

Figura 1 – Esquema de central hidrelétrica de baixa queda com turbina hélice do tipo

tubular-S (retirado de Quantz, 1976)-------------------------------------------------------9

Figura 2 – Modelo reduzido de turbina hélice do tipo tubular, com tubo de sucção

desmontado para visualização do rotor (LHPCH – UNIFEI)----------------------------9

Figura 3 – Esquema da seção meridional da turbina hélice considerada ----------------------------10

Figura 4 – Exemplo de rotor hélice (LHPCH – UNIFEI) ---------------------------------------------11

Figura 5 – Linhas médias dos perfis de pá ARC em uma estação radial com as

geometrias de projeto do rotor---------------------------------------------------------------12

Figura 6 – Geometria da pá ARC-------------------------------------------------------------------------13

Figura 7 – Geometrias de projeto em uma estação radial (grades do distribuidor e rotor)--------15

Figura 8 – Parametrização da geometria do rotor no trabalho de Lipej (2004)---------------------16

Figura 9 – Parametrizações parabólicas para a geometria do rotor de Souza (1989) --------------17

Figura 10 – Esquema de parametrização de uma grandeza genérica das pás, y --------------------17

Figura 11 – Esquema da seção meridional da turbina hélice considerada---------------------------20

Figura 12 – Convenção de pontos na seção meridional da turbina-----------------------------------21

Figura 13 – Representação das linhas de corrente absolutas instantâneas em uma

seção cilíndrica-------------------------------------------------------------------------------22

Figura 14 – Componentes de velocidade nas grades do distribuidor--------------------------------22

Figura 15 – Triângulos de velocidade nas grades do rotor --------------------------------------------22

Figura 16 – Elemento de fluido no recinto entre o distribuidor e o rotor----------------------------23

Figura 17 – Esquema iterativo para o cálculo dos perfis de velocidade após o estator------------31

Figura 18 – Esquema iterativo para o cálculo dos perfis de velocidade após o rotor--------------34

Figura 19 – Desvio angular do escoamento na saída de uma grade do distribuidor ---------------35

Figura 20 – Desvio angular do escoamento na saída de uma grade do rotor------------------------37

Figura 21 – Esquema geral da seqüência de cálculos para a análise do escoamento --------------48

Figura 22 – Fluxograma do projeto otimizado usando-se SQP-

fmincon --------------------------53

Figura 23 – Interpolações quadráticas do CRS6/CRSI ------------------------------------------------57

Figura 24 – Função bidimensional de Price (1977) ----------------------------------------------------63

Figura 25 – Fluxograma do projeto otimizado usando-se CRSA-------------------------------------64

Figura 26 – Esquema da seção meridional da turbina hélice projetada por Souza (1989) --------66

Figura 27 – Perfis de velocidade no projeto inicial ----------------------------------------------------69

Figura 28 – Distribuições de rc

no projeto inicial-----------------------------------------------------69

Figura 29 – Perdas no projeto inicial---------------------------------------------------------------------71

Figura 30 – Trabalho específico no projeto inicial -----------------------------------------------------71

Figura 31 – Ângulo de montagem (a) e arqueamento (b) nos projetos inicial e ótimo ------------75

Figura 32 – Perfis otimizados para as pás do rotor-----------------------------------------------------76

Figura 33 – Perfis de velocidade--------------------------------------------------------------------------76

Figura 34 – Perfis de momento angular por unidade de massa do fluido ---------------------------76

Figura 35 – Variação radial das perdas hidráulicas ----------------------------------------------------77

Figura 36 – Distribuição radial de trabalho específico-------------------------------------------------78

Figura 37 – Histórico de otimização da solução ‘CRSI - 3’ em termos de

(a) iterações e

(b) chamadas da função----------------------------------------------------------------------83

Lista de Tabelas

Tabela 1 − Conjunto de modelos de perdas para a turbina--------------------------------------------44

Tabela 2 − Comparação de três versões de CRSA -----------------------------------------------------62

Tabela 3 – Características básicas do projeto inicial---------------------------------------------------66

Tabela 4 – Valores adotados para os fatores empíricos de perdas------------------------------------67

Tabela 5 – Dados de projeto, resultados ótimos do ensaio e resultados do programa

para o projeto inicial de Souza (1989)-----------------------------------------------------67

Tabela 6 – Análise do projeto inicial---------------------------------------------------------------------68

Tabela 7 – Ponto de projeto e restrições de altura para a otimização --------------------------------72

Tabela 8 – Restrições de faixa para as variáveis de projeto-------------------------------------------73

Tabela 9 – Três soluções encontradas usando-se o SQP-

fmincon ----------------------------------73

Tabela 10 – Soluções usando-se o CRSI modificado --------------------------------------------------81

Tabela 11 – Comparação do projeto inicial com os melhores resultados de otimização----------83

Simbologia

Letras Latinas

a Coeficiente quadrático da parametrização

Aceleração

b Largura axial da pá, coeficiente linear da parametrização

B Altura radial da pá, aproximação para a matriz Hessiana

c Termo independente da parametrização

c Velocidade absoluta

D Diâmetro

f Função objetivo, arqueamento do perfil da pá ARC à meia corda

Força

FE Número de chamadas da função objetivo

g Função de restrição, coordenada do centróide, aceleração da gravidade

h Pior ponto da população

H Altura de queda (energia por unidade de peso)

Integral para o equilíbrio radial após o distribuidor

Integral para o equilíbrio radial após o rotor

a K

Coeficientes das regressões para rc

l Comprimento da corda do perfil da pá

l Melhor ponto da população

m Função na correlação de desvio, massa

M Fator de penalização

n Velocidade de rotação do rotor, número de variáveis de projeto

Rotação específica referente à vazão

vii

N Número de estações radiais, tamanho da população

pá

Número de pás

p Pressão, direção de busca do SQP

p Ponto tentativa

P Potência, conjunto de pontos da população

Q Vazão volumétrica

r Raio genérico, coordenada radial no plano meridional, coordenadas dos pontos

tomados no CRSA

R Coordenada radial adimensional, raio de curvatura

Número de Reynolds

ℜ Conjunto dos números reais

S Região viável do problema de otimização

t Passo da grade

TS Sucessos totais

u Velocidade circunferencial de um ponto de raio r da pá (u =

V Velocidade genérica

w Velocidade relativa

x Vetor de variáveis de projeto

x Comprimento na pá ARC

, X

Coeficientes de perda no tubo de sucção

y Grandeza genérica a ser parametrizada

Y Trabalho específico, energia por unidade de massa

z Coordenada axial

Letras Gregas

Ângulo do escoamento absoluto, ângulo geométrico e de montagem da palheta

diretriz, medida de variabilidade local

Ângulo do escoamento relativo, ângulo geométrico e de montagem da pá

Ângulo de desvio do escoamento na saída da grade

Deflexão do escoamento na grade, tolerância na população

viii

Ângulo na pá ARC

Rendimento

Coeficiente de perda de perfil, parâmetro na correlação de Soderberg

Ângulo polar, ângulo na pá ARC, coordenada circunferencial

Coeficiente de perda por choque, fator de sub-relaxação

Viscosidade absoluta

π 3,14159265

Massa específica

Ângulo de arqueamento do perfil, ângulo na pá ARC

Velocidade angular do rotor,

= 2πn

∇ Operador nabla

Superescritos

* Referente ao ponto ótimo, restrição ativada

atual Distribuição de c

atual

grade Distribuição de c

calculada de acordo com as relações de grade

L Limite inferior

n Dimensão do espaço

novo Nova distribuição de c

T Transposto

U Limite superior

Subscritos

0 Deflexão nula, ponto de partida

1 Seção de entrada do distribuidor, deflexão não-nula

2 Seção de saída do distribuidor

4 Seção de entrada do rotor

5 Seção de saída do rotor

ch Componente de choque

d Distribuidor

e Eixo

f Referente ao fluido

h Cubo ou raiz da pá, pior valor objetivo, hidráulico

inc Componente de incidência

j j-ésima coordenada do vetor

l Melhor valor objetivo

L Limite inferior

m Meio da pá, componente meridional

mec Mecânico

p Referente ao ponto tentativa

pá Pá ou rotor

P Perda

Pch Perda por choque

Pd Perda no distribuidor

Pr Perda no rotor

P(r + ts) Perda no rotor e no tubo de sucção

Pts Perda no tubo de sucção

r Componente radial, rotor

S Condições de estagnação

t Ponta da pá

u Componente circunferencial

U Limite superior

z Componente axial

Abreviaturas

ANEEL Agência Nacional de Energia Elétrica

ARC Referente à pá em formato de arco de circunferência

BFGS Broyden, Fletcher, Goldfarb e Shanno

CERPCH Centro Nacional de Referência em Pequenas Centrais Hidrelétricas

CFD Dinâmica dos fluidos computacional – Computational Fluid Dynamics

CPU Central única de processamento

CRSA Algoritmo de Busca Aleatória Controlada – Controlled Random Search Algorithm

PC Computador pessoal

SQP Programação Quadrática Seqüencial – Sequential Quadratic Programming

Siglas

LHPCH

Laboratório Hidromecânico para Pequenas Centrais Hidrelétricas

IEM

Instituto de Engenharia Mecânica

UNIFEI

Universidade Federal de Itajubá

Capítulo 1

INTRODUÇÃO & MOTIVAÇÃO

Turbinas hidráulicas são máquinas de fluxo de longa história. Vem sendo projetadas,

construídas e colocadas em operação há cerca de duzentos anos. Inicialmente como

substitutas das milenares rodas d’água no acionamento de moinhos, teares e pequenas

manufaturas, as turbinas hidráulicas são hoje, quase que exclusivamente, destinadas à geração

de energia elétrica. Para este fim específico, aliás, qualquer motor hidráulico moderno deve

preencher os seguintes requisitos técnicos básicos (Quantz, 1976): 1º) devem possibilitar o

aproveitamento de uma grande gama de saltos, cobrindo ampla faixa de alturas e vazões

disponíveis; 2º) o aproveitamento deve efetuar-se com bons valores de rendimento e com boas

características hidrodinâmicas, permitindo o acoplamento do motor hidráulico às máquinas

geradoras ainda que sejam variáveis as condições do salto (altura e vazão), de modo que a

instalação seja rentável; 3°) o eixo/árvore poderá dispor-se horizontal, inclinado ou

verticalmente, segundo o exija o acoplamento às máquinas geradoras; 4º) a velocidade

angular deve ser a mais elevada possível para que se consiga dessa forma acoplamentos

diretos ou transmissões com poucas multiplicações; 5º) devem apresentar boa regulagem a

fim de que sejam tão adequados quanto outros tipos de motores (turbinas a vapor e a gás,

motores diesel) para o serviço nas centrais elétricas; 6º) todos os elementos importantes,

especialmente os órgãos de regulagem e mancais, devem ser de fácil manutenção. As

modernas turbinas hidráulicas dos tipos Pelton, Francis, Kaplan e hélice cumprem bem todas

essas condições, superando largamente outros tipos de motores hidráulicos e competindo

economicamente com outras máquinas motoras, como os já citados motores diesel, turbinas a

vapor e a gás.

Nos primeiros projetos de turbinas hidráulicas, a experiência do próprio

engenheiro/projetista, juntamente com numerosos e dispendiosos testes com modelos tipo

tentativa-e-erro, constituíam as principais ferramentas de projeto disponíveis. A quantidade de

informação empírica era inclusive comparável aos métodos analíticos viáveis de até então.

Parte desse conhecimento empírico foi condensada em diversos diagramas e guias de projeto

(usados ainda hoje) que fornecem linhas gerais para o dimensionamento básico das turbinas

(Cordier, 1955; Quantz, 1976; Schweiger e Gregori, 1989). Outra parte desse conhecimento

ficou retida pelos próprios projetistas, sendo transmitida de “mão-em-mão”, como uma

herança, aos próximos times de engenheiros das empresas. De fato, em comparação com

outros tipos de máquinas de fluxo, como bombas, ventiladores, turbocompressores e turbinas

a gás, pode-se afirmar que são escassas as publicações técnicas referentes a projeto e

otimização de turbinas hidráulicas.

O desenvolvimento de computadores digitais na segunda metade do século XX e sua

aplicação à análise do escoamento em turbomáquinas, impulsionada primordialmente pelos

avanços no campo das turbinas a gás aeronáuticas (Denton, 1993), tornou possível o uso de

métodos complexos de simulação numérica de escoamentos para análise e projeto também de

turbinas hidráulicas. “O projeto hidrodinâmico e a construção de turbinas hidráulicas tem sido

mais uma arte do que uma ciência. Os elementos científicos tornaram-se mais numerosos com

os recentes avanços na tecnologia de análise de escoamento” (Ueda, 1982). Atualmente,

programas dos tipos Euler 3D e Navier-Stokes 3D já são ferramentas-padrão no

desenvolvimento de novas unidades de turbinas hidráulicas, podendo em certos casos ser até

usados com rotinas de otimização (Lipej, 2004; Peng et al., 2002a e 2002b). Detalhes da

separação do escoamento, fontes de perdas e suas distribuições em componentes, análise

acoplada de componentes no ponto de projeto e fora dele, e baixos níveis de pressão com

risco de cavitação agora são problemas mais amenos de se analisar com a assim denominada

Dinâmica dos Fluidos Computacional – CFD (Drtina e Sallaberger, 1999; “Design by

Numbers”, 1998).

A aplicação dessas técnicas modernas de CFD para a predição do campo de

escoamento através de uma turbina inteira tem levado a uma melhor compreensão física dos

fenômenos que ocorrem nesses escoamentos, com conseqüências diretas sobre o projeto

hidrodinâmico dos componentes da turbina. Além disso, o progresso nas técnicas

experimentais de medição e testes com modelos é outro fator importante que tem contribuído

para essa compreensão mais detalhada dos fenômenos fluido-dinâmicos em turbinas

hidráulicas. No tocante à parte experimental, inclusive, os avanços na análise numérica

computacional e a tecnologia de predição das características de funcionamento não eliminam

os ensaios com modelos como meio para se melhorar o rendimento, especialmente fora do

ponto de projeto. Tais ensaios, no entanto, agora podem ser muito mais objetivos, sendo

realizados em menor número e já na fase final de projeto/prototipagem, reduzindo-se

significativamente, assim, os custos com experimentos (Ueda, 1982; “Design by Numbers”,

1998).

Mas embora os programas do tipo Navier-Stokes 3D apresentem resultados bastante

confiáveis, revelando detalhes importantes do escoamento local e fornecendo predições de

performance muito precisas, um esforço computacional significativo é exigido com a geração

e modificação de malhas e com a resolução das complicadas equações de movimento

(viscoso, 3D e não-permanente) em cada investigação numérica. Por exemplo, na otimização

de turbomáquinas, o contexto desta dissertação, quando uma alteração geométrica é efetuada

pelo algoritmo de otimização, as malhas devem ser recalculadas e o campo de escoamento

deve ser novamente avaliado pelo solver – com seu alto custo computacional. Esse esforço

freqüentemente impede a integração de simulações sofisticadas do tipo Navier-Stokes 3D ao

longo de todas as etapas de projeto (Drtina e Sallaberger, 1999; “Design by Numbers”, 1998).

Além disso, uma análise de precisão razoável, mas simples e rápida, ainda é essencial para as

fases iniciais do projeto, quando a geometria não está complemente determinada (Oh e Kim,

2001; Yoon, et al., 1998). Em turbinas a gás, por exemplo, são bastante comuns as

publicações sobre métodos computacionais de baixo custo para análise e projeto preliminares.

Nesse âmbito, é típica a aplicação de Métodos de Curvatura de Linha de Corrente com uma

modelagem simplificada para as perdas e desvios do escoamento (Yoon et al., 1998; Lee e

Chung, 1991; Park e Chung, 1992; Sullerey e Kumar, 1984) ou mesmo análises apenas na

linha média – 1D – (Kacker e Okapuu, 1982; Souza Júnior et al., 2005), talvez ainda

indispensáveis para a otimização inicial de um novo projeto e para a predição dos

rendimentos atingíveis.

Tendo essas considerações em mente, pode-se afirmar que metodologias

intermediárias de projeto, com baixo custo computacional, são necessárias para as etapas

iniciais de projeto otimizado de qualquer turbomáquina. E no caso particular de turbinas

hidráulicas, são ainda raramente encontradas na literatura. Essas metodologias consistem no

que se chama de otimização conceptual, fornecendo uma geometria simples mas

suficientemente representativa para a descrição dos principais parâmetros de projeto de

rotores e estatores e também levando a tendências corretas para o campo de escoamento

ótimo. Neste trabalho, uma tal metodologia é desenvolvida para turbinas hidráulicas axiais, do

tipo hélice-tubular.

As configurações dos tipos bulbo e tubular, inclusive, vêm sendo usadas cada vez mais

em diversos países, entre os quais o Brasil, em detrimento às turbinas Kaplan convencionais

de eixo vertical. Há de fato uma nítida tendência mundial em direção aos aproveitamentos de

baixas e baixíssimas quedas (inferiores a 15 m) até então inexplorados por questões

econômicas, mas que agora, em virtude do esgotamento dos aproveitamentos tradicionais,

com quedas moderadas e altas, e por restrições ambientais cada vez mais fortes, despontam

como excelente alternativa para a expansão da matriz hidroelétrica mundial, especialmente na

América do Sul, Índia, China, Sudeste Asiático e Oeste e Sudoeste Africanos (Dansie, 1996;

Subrahmanyam, 1982, Hindley, 1996). No Brasil, em particular, o potencial hidrelétrico situa-

se ao redor de 260 GW, enquanto que a potência instalada em nossas hidrelétricas é de cerca

de 80 GW, ou seja, menos de um terço do potencial (sites do CERPCH e da ANEEL). Há

portanto muito espaço ainda para o desenvolvimento hidrelétrico brasileiro. No tocante às

turbinas axiais dos tipos bulbo e tubular, especialmente destinadas às baixas quedas, o

destaque vai para a região Amazônica, onde diversos rios de planície podem ser aproveitados.

Por exemplo, já se encontra em fase de projeto básico as duas maiores centrais de baixa queda

do mundo, a saber, as usinas hidrelétricas de Jirau e Santo Antônio no rio Madeira. São nada

menos que oitenta e oito turbinas do tipo bulbo, com rotores de 8 m de diâmetro, que

produzirão um total de 6450 MW, um valor significativo para o desenvolvimento da Região

Norte do Brasil (Porto et al., 2005).

Na metodologia proposta neste trabalho, dado o ponto de projeto (vazão e rotação),

restrições para a altura de queda e algumas dimensões preestabelecidas para a turbina, busca-

se uma configuração ótima para as geometrias do distribuidor e do rotor de modo que a

energia do escoamento seja absorvida da maneira mais eficiente, ou seja, com máximo

rendimento.

O projeto otimizado de outros tipos de máquinas de fluxo axiais, como bombas,

ventiladores, turbocompressores, turbinas a vapor e turbinas a gás também poderia ser

realizado conforme a metodologia deste trabalho, feitas as devidas modificações. A

parametrização da geometria, a condição de equilíbrio radial e as técnicas de otimização

seriam essencialmente iguais às apresentadas nos capítulos subseqüentes.

Capítulo 2

FORMULAÇÃO DO PROBLEMA DE OTIMIZAÇÃO

Neste capítulo, define-se o problema de otimização de projeto da turbina. Em linhas

gerais, consiste na busca de alguns parâmetros geométricos básicos para as pás do rotor e para

as aletas do distribuidor – as variáveis de projeto – de forma a maximizar o rendimento da

turbina – a função objetivo –, dadas algumas dimensões preestabelecidas, a vazão volumétrica

e a rotação da máquina. A altura de queda disponível, resultante da geometria de projeto

investigada e do par vazão-rotação definido, deve permanecer entre limites inferior e superior,

sendo essas as restrições não-lineares do problema. Há também restrições laterais (ou de

faixa) para as variáveis de projeto, definindo-se a região viável, ou região de busca do

problema.

Em termos mais formais, tem-se um problema de minimização não-linear, com

restrições, na forma:

minimizar f(x)

sujeito a g

(x) ≤ 0, i = 1, ..., m (2.1)

x ∈ S

x é o vetor n-dimensional das variáveis de projeto x

, j = 1, ..., n. A região de busca S é

definida por limites inferior e superior, x

e x

respectivamente, para cada componente de x:

S = {x ∈ ℜ

: x

≤ x

, j = 1, ..., n}. A função objetivo é f(x) = −

(x), em que

é o

rendimento da turbina (definida geometricamente pelo vetor x e operando no ponto vazão-

rotação dado). g

(x), i = 1, ..., m, são as m funções de restrição, a saber:

(x) = H

– H(x) (2.2)

(x) = H(x) – H

(2.3)

em que H(x) é a altura de queda disponível da turbina e H

e H

são respectivamente limites

inferior e superior, de modo que se tenha H

≤ H ≤ H

Essa forma de impor as restrições quanto à altura de queda disponível, H, é adequada

quando se usam métodos de otimização baseados em gradiente e outras derivadas direcionais,

como Programação Quadrática Seqüencial, Máximo Descenso, Gradiente Conjugado, etc.

Uma outra maneira de impor as restrições quanto à altura é penalizar a função objetivo:

⎪

⎩

⎪

⎨

⎧

>−+−

≤≤−

<−+−

HHHHM

HHH

HHHHM

,)(

(2.4)

O fator de penalização M é um número positivo suficientemente grande. Novamente, o

objetivo é maximizar

(minimizar −

) com H pertencendo ao intervalo [H

; H

]. A escolha

do fator de penalização M não deve conduzir o processo de otimização em direção a uma

minimização da penalidade apenas, perdendo-se a informação importante da função objetivo,

i.e., o rendimento

. Por outro lado, as restrições não devem ser violadas ao fim do processo

de convergência. Logo, devem-se perfazer alguns testes a fim de se poder escolher valores

satisfatórios para M.

O método de penalização é mais apropriado quando se usam métodos de otimização

dos tipos populacionais, como os Algoritmos de Busca Aleatória Controlada, Genéticos e de

Evolução Diferencial. Porém, à medida que cresce o número de restrições a serem impostas,

um esquema de penalização pode causar sérias instabilidades no processo de convergência,

chegando até mesmo a fazer fracassar a otimização. A escolha dos diversos fatores de

penalização nesses casos também torna-se complicada. O mesmo pode ocorrer em problemas

multi-objetivos tratados como mono-objetivos por meio de somas ponderadas das diversas

funções objetivo originais. Em tais casos, um outro método deve ser empregado para tratar

eficientemente as restrições e/ou os vários objetivos (Oyama et al., 2005).

Como se percebe, será feita a maximização do rendimento apenas no ponto de projeto,

caracterizado pelo par vazão-rotação dado e pela faixa de alturas de queda imposta como

restrição. No entanto, como regra geral, pode-se afirmar que, ao se aumentar o rendimento no

ponto de projeto (comumente, o rendimento máximo), melhora-se a eficiência da máquina em

uma ampla faixa operacional em torno desse ponto (Ueda, 1982). Assim, espera-se que as

turbinas otimizadas segundo a metodologia desenvolvida neste trabalho apresentem boas

características de funcionamento também fora do ponto de projeto, mesmo que isso não tenha

sido colocado como objetivo do problema de otimização.

Outro ponto que pode (e deve) ser levantado é o fato de não se inserir algum

parâmetro de cavitação nos objetivos ou nas restrições do problema. De fato, embora a

avaliação do fenômeno de cavitação seja um aspecto básico no projeto de turbinas hidráulicas

(Lipej, 2004; Peng et al., 2002a), este trabalho está focado apenas no rendimento da turbina.

Isso, no entanto, ainda é viável para um projeto realista porque a ocorrência de cavitação pode

ser evitada preliminarmente controlando-se, por exemplo, os ângulos de incidência do

escoamento na entrada do rotor e o carregamento hidrodinâmico sobre as pás. Especificando-

se critérios como esses, pode-se realmente proceder à otimização do rendimento da turbina

com alguma segurança contra os riscos de cavitação, o que é aceitável num estágio inicial de

projeto. Como complemento a essa metodologia intermediária, uma técnica mais sofisticada,

usando-se CFD, pode ser aplicada no contexto de um projeto inverso otimizado para as grades

do distribuidor e do rotor. Nessa etapa mais avançada, os perfis de velocidade inicialmente

determinados na otimização conceptual são preservados ou levemente alterados ao mesmo

tempo em que os perfis das pás são desenhados de modo a garantir mínima ocorrência de

cavitação.

Capítulo 3

PARAMETRIZAÇÃO DA GEOMETRIA DAS PÁS

3.1 MODELO GEOMÉTRICO

A turbina hidráulica axial considerada neste trabalho é uma turbina hélice do tipo

tubular, como mostram as Figuras 1, 2 e 3. Admite-se ainda um distribuidor cilíndrico

(i.e., sem conicidade) e aletas sem torção ao longo de sua envergadura.

O tratamento da geometria do rotor pode ser feito de diferentes formas,

dependendo das informações que se pretendem obter a partir dessa geometria. Por

exemplo, para um estudo detalhado do carregamento hidrodinâmico sobre as pás, com

identificação de pontos de pressão mínima e pontos de descolamento da camada-limite,

seria necessária uma descrição pormenorizada dos perfis das pás em qualquer estação

radial. Um método de painéis, por exemplo, talvez fosse satisfatório para essa descrição

em algumas seções cilíndricas, e a geometria assim discretizada poderia ser interpolada

de alguma maneira para que a descrição fosse completa ao longo de toda a extensão

radial das pás (

span

A quantidade de parâmetros e variáveis de projeto escolhidos deve ser

compatível como o nível de informação que se pretende gerar a partir da geometria da

turbina. Um aspecto fundamental é a seleção dessas variáveis, que deve ser feita de

modo a facilmente se identificarem as melhorias de desempenho na fase inicial de

projeto. O número de variáveis de projeto também merece atenção, pois influencia

significativamente a convergência do processo de otimização.

Figura 1 – Esquema de central hidrelétrica de baixa queda com turbina hélice do tipo

tubular-S (retirado de Quantz, 1976).

Figura 2 – Modelo reduzido de turbina hélice do tipo tubular, com tubo de sucção

desmontado para visualização do rotor (LHPCH – UNIFEI).

Para o nível de otimização a que se propõe este trabalho, será suficiente a

avaliação das distribuições de velocidade apenas nas seções de entrada e de saída do

distribuidor e do rotor. Esse tipo de informação já viabiliza a estimativa dos principais

parâmetros energéticos da turbina, como potências, perdas e rendimentos. Não será

necessário, por exemplo, o cálculo do campo de escoamento ao redor dos perfis das pás,

posto que a principal preocupação é com relação às deflexões resultantes do

escoamento nas grades. Logo, a geometria será simplificada de forma que, juntamente

com as equações da continuidade, equilíbrio radial e correlação de desvio, seja a

estritamente necessária para o cálculo das distribuições de velocidade nas seções de

entrada e de saída do distribuidor e do rotor.

Tubo de Sucção

Rotor Distribuidor

Eixo

Figura 3 – Esquema da seção meridional da turbina hélice considerada.

Turbinas hidráulicas axiais, especialmente as destinadas às baixas quedas, como

as dos tipos bulbo e tubular, apresentam pás muito pouco arqueadas. A Figura 4 é a

fotografia de um rotor hélice para alturas de queda nem tão baixas, visto que ele

apresenta cinco pás; os rotores de baixas e baixíssimas quedas tem apenas quatro ou até

mesmo três pás. Mesmo assim, acham-se na Figura 4 perfis delgados e pouco arqueados

para as pás do rotor, com exceção talvez para as estações mais próximas ao cubo (que

são mais espessas devido a requisitos estruturais).

Decidiu-se então aproximar as linhas médias dos perfis das pás (

camber lines

)

por arcos de circunferência (ARC) de pequena curvatura, uma escolha razoável para um

rotor de turbina hidráulica axial de baixa queda. A espessura das pás não será

considerada neste trabalho, pois, como já se justificou, nenhum fenômeno de cavitação

ou separação do escoamento será avaliado diretamente pela metodologia adotada. De

fato, quando os perfis são suficientemente delgados, a espessura não contribui para a

força de sustentação sobre o hidrofólio: em um perfil delgado e pouco arqueado, a

espessura afeta apenas a distribuição de pressão, a sustentação sendo uma função

somente do ângulo de ataque e do arqueamento do perfil (teoria do aerofólio delgado –

Karamcheti, 1980). Será, pois, suficiente considerar apenas as linhas médias (ou de

arqueamento) para o cálculo das velocidades através da turbina, sem considerar as

espessuras.

Figura 4 – Exemplo de rotor hélice (LHPCH – UNIFEI). Notar a esbeltez dos perfis.

Dessa forma, em uma dada estação radial, o perfil da pá ARC é completamente

definido pelo ângulo de montagem da corda (

stagger angle

, comprimento da corda,

, e arqueamento máximo (posicionado à meia corda),

, Figura 5. Como o número de

pás é previamente fixado neste trabalho, ou seja, não é uma variável de projeto, a

especificação da corda e da flecha do perfil é equivalentemente obtida pela

especificação das razões corda-passo (

chord-pitch ratio

, e flecha-corda (

relative

camber

, consideradas mais apropriadas de se trabalhar por serem adimensionais.

Portanto, os parâmetros geométricos de projeto de uma grade do rotor são,

Essa escolha é adequada porque, além de serem grandezas diretas para se especificar

num dimensionamento, levam a todas as características geométricas e cinemáticas

necessárias das grades numa fase inicial de projeto, tais como ângulos de incidência,

ângulos de arqueamento, ângulos de ataque, ângulos de desvio, deflexões, etc.

Figura 5 – Linhas médias dos perfis de pá ARC em uma estação radial com as

geometrias de projeto do rotor.

A seguir, mostra-se como, a partir de

, obtêm-se os demais parâmetros

geométricos da pá ARC. Numa dada estação radial de raio

, o passo

é dado por:

pá

(3.1)

onde

pá

é o número de pás. Logo, dadas as relações

, calculam-se os

comprimentos da corda e da flecha respectivamente por:

tt)/(ll =

(3.2)

ll)/( ff =

(3.3)

definindo-se completamente a linha média da pá ARC, Figura 5. A corda axial

é dada

por:

senl=b

(3.4)

Tendo em vista a Figura 6, tem-se pelo teorema de Pitágoras:

222

)2/(

pá

Rx =+l

(3.5)

A flecha

é tal que:

xRf

pá

−=

(3.6)

pelo que:

)( fRx

pá

−=

(3.7)

Das Equações (3.5) e (3.7), obtém-se:

2222

24/

pápápá

RffRR =+−+l

(3.8)

pá

(raio de curvatura do perfil) (3.9)

equação essa que mostra ser o raio de curvatura do perfil,

pá

, função apenas da corda

e da flecha

(

)

Figura 6 – Geometria da pá ARC.

Tem-se ainda:

⎟

⎠

⎞

⎜

⎝

⎛

pá

arcsen

(3.10)

−= 2/

(3.11)

−= 2/

(3.12)

Com isso,

−

−=

(ângulo geométrico de entrada da pá) (3.13)

Para uma pá em ARC, o raio de curvatura também pode ser calculado por:

coscos

ββ

−

pá

(3.14)

de modo que

⎟

⎠

⎞

⎜

⎝

⎛

cosarccos

ββ

pá

(ângulo geométrico de saída da pá) (3.15)

As Equações (3.1) a (3.4), (3.9), (3.13) e (3.15) fornecem todas as grandezas

geométricas necessárias para o cálculo do escoamento nas seções de entrada e de saída

das grades do rotor.

Para a geometria do distribuidor, como as aletas não são torcidas ao longo do

raio, um único ângulo de saída,

, é suficiente como variável de projeto, devendo o

restante da geometria das aletas ser previamente estabelecido. Essa escolha permitirá ao

programa de otimização alterar a quantidade de movimento angular do fluido que o

distribuidor “descarrega” para o rotor, buscando-se um valor tal que favoreça as

melhores condições de absorção da energia do fluido pelas pás. De fato, a otimização

da geometria somente do rotor, em geral, é insuficiente para se obter grandes aumentos

no rendimento da turbina como um todo (Ueda, 1982).

Na Figura 7, apresentam-se as geometrias de projeto numa estação radial para as

grades do distribuidor e do rotor.

Figura 7 – Geometrias de projeto em uma estação radial (grades do distribuidor e rotor).

3.2 PARAMETRIZAÇÃO

Escolheu-se parametrizar o ângulo de montagem,

, a razão corda-passo,

, e o

arqueamento relativo,

, em termos dos valores dessas grandezas nas estações radiais

referentes ao cubo, meio e ponta das pás, ou seja, 0%, 50% e 100% da envergadura,

respectivamente. Isso leva a 3

3 = 9 variáveis de projeto para o rotor. Somando ainda

a única variável de projeto escolhida para o distribuidor (

), resulta um total de apenas

9 + 1 = 10 variáveis de projeto para a turbina inteira.

Em geral, uma parametrização deve permitir boa reprodução de toda a

geometria, ainda que alguns detalhes de menor importância sejam negligenciados. Além

de facilitar as alterações geométricas durante o processo de otimização, a

parametrização é muito importante para se reduzir o número de variáveis de projeto

numa turbomáquina, sendo usada até mesmo na geração e modificação de malhas em

problemas tratados via CFD (Drtina e Sallaberger, 1999; Lipej, 2004; Oyama, 2005).

Dentre as muitas possibilidades de parametrização, optou-se por funções

parabólicas para se descrever a geometria das pás ao longo de toda sua extensão radial.

Há de fato diversas maneiras de se parametrizar a geometria, sendo que a forma mais

adequada e conveniente talvez só possa ser determinada por meio de um estudo

paramétrico. Neste trabalho, tal investigação não foi realizada, mas mesmo assim

acredita-se que a escolha de parábolas seja razoável para aproximar satisfatoriamente as

configurações geométricas usualmente encontradas em turbinas hidráulicas axiais. Na

Figura 8, reproduzida do artigo de Lipej (2004), encontra-se o ângulo de montagem (

)

e o arqueamento relativo (

) segundo as parametrizações adotadas por esse autor.

Embora não as cite explicitamente, dizendo somente que poderiam ser compostas por

funções lineares ou de ordem mais alta, percebe-se claramente pela figura que funções

polinomiais, exponenciais ou logarítmicas são fortes candidatas a servir de

parametrização.

Figura 8 – Parametrização da geometria do rotor no trabalho de Lipej (2004).

Além disso, a geometria do rotor da turbina hidráulica projetada por Souza

(1989) – que será usada mais adiante para as comparações com os resultados da

otimização – é razoavelmente reproduzida pelas parametrizações parabólicas propostas,

como mostra a Figura 9.

Seja

qualquer uma das três grandezas a serem parametrizadas (

são respectivamente os valores de

nas estações do cubo, meio e ponta das pás,

cujos raios, adimensionalizados pelo raio da ponta (

), são

) e

= 1), Figura 10. Tem-se então o seguinte sistema linear em

oriundo

da parametrização

² +

≤

1 :

⎪

⎩

⎪

⎨

⎧

=++

mmm

hhh

ycba

ycbRaR

(3.16)

0.0 0.2 0.4 0.6 0.8 1.0

Raio adimensional, r/rt

Ângulo de montagem da corda, graus

Projeto de Souza (1989)

Parametrização

0.00.20.40.60.81.0

Raio adimensional, r/ rt

0.0

0.2

0.4

0.6

0.8

1.0

1.2

1.4

1.6

1.8

Razão corda-passo, (-)

Projeto de Souza (1989)

Parametrização

0.0

0.2

0.4

0.6

0.8

1.0

1.2

1.4

1.6

1.8

Figura 9 – Parametrizações parabólicas para a geometria do rotor de Souza (1989).

= 1

Figura 10 – Esquema de parametrização de uma grandeza genérica das pás,

A forma matricial para o sistema (3.16) é:

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

111

( [

]{

abc

} = {

} ) (3.17)

Como

(=1) são todos distintos entre si, a matriz [

] é não-singular pelo

teorema de Vandermonde e logo esse sistema pode ser resolvido pela regra de Cramer

(Alonso

et al

., 1990). Para tal, calcula-se:

)1)((

111

det

−−+−==

hmhmhmmm

RRRRRRRR

(3.18)

)()1()1(

det

mhthmmh

RRyRyRy

a −+−+−==

(3.19)

)()1()1(

det

22222

hmthmmh

RRyRyRy

b −+−+−==

(3.20)

Agora, pela regra de Cramer:

det

(3.21)

det

(3.22)

e, pela terceira equação do sistema (3.16):

bayc

−

−=

(3.23)

o que conclui o cálculo da parábola de parametrização adotada.

Não é exagero enfatizar, novamente, o caráter conceptual do projeto otimizado

proposto neste trabalho. Nesse estágio do dimensionamento da turbina, os padrões

obtidos para os perfis de velocidade a partir da geometria considerada são mais

relevantes que a própria geometria, embora esta também seja representativa

construtivamente. Um procedimento complementar poderia ser a otimização do projeto

inverso das grades, garantindo as distribuições de velocidade inicialmente calculadas na

otimização conceptual. De fato, tal abordagem vem sendo efetuada em questões mais

complexas de otimização de turbinas hidráulicas, onde o desempenho de cavitação

também é um objetivo do problema e programas CFD do tipo Navier-Stokes 3D são

empregados como

solver

. No entanto, como já citado anteriormente, o alto custo

computacional nesse tipo de análise ainda é um fator limitante para o processo de

otimização. Nesta dissertação, mais uma vez, o objetivo é desenvolver uma

metodologia de projeto intermediária, de baixo custo computacional, viável de se

perfazer rapidamente em um único computador (PC) e que preceda a análises mais

sofisticadas para já adiantar resultados básicos importantes.

Capítulo 4

CÁLCULO DO ESCOAMENTO EM UMA TURBINA

HIDRÁULICA AXIAL

Repete-se a seguir a Figura 3, onde se apresenta o esquema da seção meridional

da turbina hidráulica considerada neste trabalho, i. e., uma turbina hélice do tipo

tubular, com distribuidor cilíndrico e aletas sem torção ao longo do raio.

Figura 11 – Esquema da seção meridional da turbina hélice considerada.

Distribuidor Rotor

Tubo de sucção

Neste capítulo, será desenvolvido o procedimento de cálculo das distribuições de

velocidade antes e após o distribuidor (estator) e antes e após o rotor. Para tal, será

imposta a condição de equilíbrio radial das máquinas de fluxo axiais. O desvio angular

do escoamento na saída das grades será avaliado usando-se correlações empíricas

provenientes da literatura. Uma vez estabelecidas as distribuições de velocidade,

calculam-se as perdas hidráulicas através da turbina, também por meio de correlações, e

a potência absorvida pelas pás, usando-se agora a equação de Euler para as

turbomáquinas – equação integral da quantidade de movimento angular. Com isso, tem-

se o rendimento da turbina (o alvo da otimização) e a altura de queda disponível, que

será submetida a restrições.

Os desenvolvimentos que se seguem fazem referência às Figuras 12 a 15 para a

convenção de pontos e definições dos triângulos de velocidade. Os índices

referem-se respectivamente às estações radiais do cubo (

hub

) e da ponta das pás (

tip

é a velocidade absoluta,

é a velocidade relativa e

é a velocidade circunferencial (ou

de condução) das pás.

são respectivamente os ângulos do escoamento absoluto e

relativo, medidos a partir da direção circunferencial (convenção alemã). Os índices

correspondem respectivamente aos componentes circunferencial e meridional da

velocidade. Por exemplo,

é o componente meridional da velocidade absoluta na

seção de saída do rotor, na estação do cubo.

Figura 12 – Convenção de pontos na seção meridional da turbina.

1: entrada do distribuidor; 2: saída do distribuidor

4: entrada do rotor; 5: saída do rotor

Fluxo Tubo de sucção

5 4 2 1

Distribuidor Rotor

Figura 13 – Representação das linhas de corrente absolutas instantâneas

em uma seção cilíndrica.

Figura 14 – Componentes de velocidade Figura 15 – Triângulos de velocidade

nas grades do distribuidor. nas grades do rotor.

RotorDistribuidor

4.1 EQUAÇÃO DE EQUILÍBRIO RADIAL

A equação da quantidade de movimento na direção radial é também conhecida

como equação de equilíbrio radial. Em máquinas de fluxo axiais, essa equação

desempenha papel fundamental no estabelecimento das distribuições de velocidade.

Neste trabalho, a equação de equilíbrio radial será usada em conjunto com as equações

da energia e da continuidade para relacionar os componentes meridional e

circunferencial da velocidade absoluta após o distribuidor (

) e após o rotor (

). Isso é realmente necessário para que o campo de velocidade a ser obtido seja

compatível com os princípios físicos de conservação de energia, massa e quantidade de

movimento.

A equação de equilíbrio radial pode ser obtida da maneira como se segue.

Admite-se um escoamento com simetria axial e incompressível. Desprezam-se as forças

viscosas e o componente radial da velocidade, ou seja, supõe-se um fluxo puramente

axial na projeção meridional. Admite-se ainda escoamento absoluto em regime

permanente no estator e escoamento relativo estacionário no rotor.

Figura 16 – Elemento de fluido no recinto entre o distribuidor e o rotor.

Aplicando-se a 2ª lei de Newton para o movimento ao elemento de fluido

representado na Figura 16, tem-se:

dmd aF

(4.1)

dθ

p + dp

Distribuidor Rotor

sendo

a força resultante sobre o elemento de massa

, sua aceleração medida

de um referencial inercial. Tomando o componente radial de ambos os vetores:

)

dmd

(4.2)

O componente radial da aceleração é a aceleração centrípeta devida à velocidade

circunferencial

. Logo:

)

F −=

(4.3)

O componente radial da força resultante tem origem nas pressões estáticas

nas faces do elemento (Figura 16):

)

dLdrdpdLdrdppdLdrpd

θθθ

−=+−= )(F

(4.4)

Das Equações (4.3) e (4.4), e sendo

, tem-se:

(Equação de equilíbrio radial) (4.5)

Essa é a equação da quantidade de movimento na direção

, chamada de equação de

equilíbrio radial.

Usando agora a equação de Bernoulli para o rotor (Bran e Souza, 1979):

constante

=−+

uwp

(4.6)

onde o efeito gravitacional já fôra negligenciado, tem-se por diferenciação:

(4.7)

Observe que, como se desprezou anteriormente o componente radial da velocidade,

admitiu-se implicitamente linhas de corrente paralelas ao eixo da máquina (razoável

como primeira aproximação numa turbina hidráulica axial), e, portanto,

suposto constante ao longo de uma dada linha de corrente.

Dos triângulos de velocidade e da equação de Euler das turbomáquinas, tem-se:

)(

wwcc

ccuY

uupá

−

=−=

(4.8)

)(

dwdc

ucddY

upá

+−=−=

(4.9)

Impondo-se agora

pá

= constante na Equação (4.9)

−

trabalho específico uniforme ao

longo do raio

−

, obtém-se:

pá

(4.10)

constante=

(“vórtice potencial”) (4.11)

dwdc =

(4.12)

Além disso,

uuuuu

rdcdrcrdcdrcrcd

−

⇒=+= 0)(

(4.13)

Combinando as equações de equilíbrio radial (4.5) e Bernoulli (4.7):

drc

−=

(4.14)

e usando as Equações (4.13) e (4.12), chega-se a:

)(

rdc

−=−

(4.15)

ou seja:

0)(

=−

ccd

(4.16)

constante0

=⇔=

cdc

(4.17)

Assim,

pá

= const.

⇒

= const. e

= const. numa máquina de fluxo

hidráulica puramente axial.

Inversamente, supondo

= const., pode-se provar que

pá

= const. Da Equação

(4.9):

ducudcucddY

uuupá

−

−=−= )(

(4.18)

drcrdcdY

uupá

−−=

(4.19)

pelo que:

pá

rdc

drc −−=

(4.20)

em que

, com

sendo a velocidade angular do rotor.

Usando novamente a Equação (4.14) (equilíbrio radial e Bernoulli) e a Equação

(4.9), deriva-se:

páu

drc

−−=

(4.21)

que com a Equação (4.20) torna-se:

páu

pá

rdc

−−=

⎟

⎠

⎞

⎜

⎝

⎛

−−

(4.22)

pá

dY +−=

⎟

⎠

⎞

⎜

⎝

⎛

−

(4.23)

Usando-se agora a hipótese de

= const. na Equação (4.23), hipótese esse que leva

diretamente a

−

² = const., ou

−

² = 0, obtém-se:

01 =

⎟

⎠

⎞

⎜

⎝

⎛

−

pá

(4.24)

Logo, ou

pá

= const. O caso

é um caso particular sem interesse

prático para turbinas hidráulicas axiais, pois leva a um grau de reação teórico de apenas

0,5 e um ângulo de entrada no rotor

= 90°, valor muito alto para uma turbina axial.

Conclui-se dessa maneira que

= const.

⇔

pá

= const. (

= const.) em uma

máquina de fluxo hidráulica axial, de acordo com a condição de equilíbrio radial. A

hipótese de vórtice potencial é de fato muito usada nos projetos de turbinas axiais, pois

constitui uma solução extremamente simples para a equação de equilíbrio radial, com

distribuições uniforme para a velocidade meridional e trabalho específico, e hiperbólica

para a velocidade circunferencial. No entanto, só é verificada na prática se a torção

radial das palhetas diretrizes (aletas do distribuidor) realmente assegurar a formação do

“vórtice potencial” no recinto entre o distribuidor e o rotor. Para isso, o ângulo de saída

deve crescer do cubo para a ponta de maneira tal que, juntamente com a velocidade

meridional e o desvio angular do escoamento na saída das grades em cada estação

radial, se obtenha

= const. Se esse não for o caso, a distribuição de velocidade

meridional não será uniforme ao longo do raio, assim como a distribuição do trabalho

específico das pás. Além disso, a condição de vórtice potencial não constitui,

necessariamente, a configuração ótima para o campo de escoamento. Logo, o correto

equilíbrio radial deverá ser considerado no intuito de se obter perfis de velocidade mais

realistas e otimizados.

4.2 EQUILÍBRIO RADIAL APÓS O DISTRIBUIDOR

Admitindo escoamento incompressível e em regime permanente no distribuidor

(estator), a equação da energia (1ª lei da Termodinâmica), é escrita como (Fox e

McDonald, 1998):

PdSS

Ypp

=−

(4.25)

em que

é a pressão de estagnação (

é a massa específica do fluido e

é a perda de energia mecânica por unidade de massa que ocorre no distribuidor,

devido às irreversibilidades do escoamento.

Na entrada do distribuidor, a pressão estática,

, e a velocidade,

, são

admitidas constantes, de modo que

é uniforme ao longo do raio. Portanto, derivando

a Equação (4.25), obtém-se:

PdS

=−

(4.26)

A equação de equilíbrio radial (4.5) na saída do distribuidor é escrita como:

(4.27)

Da definição de pressão de estagnação, pode-se escrever:

)(

222 umS

ccpp ++=

(4.28)

S 2

ρρ

++=

(4.29)

Usando a equação da energia na forma (4.26) e a equação de equilíbrio radial (4.27) em

(4.29), obtém-se:

uPd 2

ρρρρ

++=−

(4.30)

que pode ser rearranjada como:

rcc

Pdm

+++=

(4.31)

rcd

Pdm

++=

)(

(4.32)

A Equação (4.32) expressa a condição de equilíbrio radial após o estator.

Distribuições de

que satisfazem a essa equação estão em conformidade com os

princípios de conservação de energia e quantidade de movimento. Embora algumas

simplificações tenham sido feitas para se chegar à Equação (4.32) (forças viscosas

desprezíveis na direção radial, fluxo puramente axial na projeção meridional, pressão de

estagnação uniforme na entrada do distribuidor), a observância dessa relação conduz a

tendências realistas para o núcleo do escoamento em turbomáquinas axiais. Novamente

pode-se observar que com a hipótese de vórtice potencial (

= const.) e desprezando-

se o efeito das perdas sobre o equilíbrio radial (

/ dr

= 0), obtém-se de (4.32)

const., tal qual fora demostrado na seção anterior usando-se a equação de Bernoulli

(que de fato pode-se aplicar quando são desprezados os efeitos viscosos).

Dividindo agora ambos os membros da Equação (4.32) por

e reagrupando,

chega-se finalmente a:

rcd

Ycd

uuPdm

)()2/(

−=

(4.33)

rcd

Ycd

uPdm

)(

)2/(

−=

(4.34)

equação essa numa forma conveniente para ser integrada e expressar

em termos de

. Como última simplificação, o efeito das perdas hidráulicas sobre o equilíbrio

radial pode ser negligenciado na Equação (4.34), ou seja, admite-se

/ dr

= 0. Note

que não está sendo afirmado serem nulas as perdas no distribuidor, mas apenas que seu

efeito sobre o equilíbrio radial na Equação (4.34) é desprezível (Peng

et al

., 2002a).

Essa afirmação será tanto mais válida quanto mais uniforme for o perfil de

. Além

disso, citando Denton (1993), “ ...muitos escoamentos são dominados pelas leis de

conservação de massa, energia e quantidade de movimento, e não por efeitos viscosos

detalhados. O poder das equações de conservação nunca deve ser subestimado e

modelos de escoamento que não satisfaçam a essas equações estão fadados ao

fracasso”. Dessa forma, tem-se:

rcd

)(

−=

(4.35)

que após a integração do cubo (

) a um raio

genérico torna-se:

∫

−=−

hmm

rcd

crc

)(

(4.36)

Definindo:

∫

−=

rcd

)(

(4.37)

escreve-se explicitamente a distribuição de velocidade meridional após o distribuidor

satisfazendo ao equilíbrio radial ao longo de todo o raio:

)()(

rIcrc

dhmm

(4.38)

A continuidade geral é agora imposta para o cálculo da velocidade meridional na

estação do cubo (

);

é a vazão volumétrica através da turbina:

∫

rdrrcQ )(2

(4.39)

)(

rdrrc

∫

(4.40)

Substituindo a Equação (4.38), obtém-se finalmente:

)(

rdrrIc

dhm

=⋅+

∫

(4.41)

A Equação (4.41) é uma equação não-linear – com a incógnita em uma integral – para

ajustar

de acordo com a continuidade geral e a condição de equilíbrio radial. Para

resolvê-la, escolheu-se um método de bisseção padrão do MatLab

(função

fzero

). Já

a integração em (4.41) é realizada numericamente segundo a regra de Simpson (Silva e

Miyazima, 2000).

O cálculo da integral

(

) em (4.37) exige o conhecimento prévio da distribuição

do momento cinemático

(

). Porém, os componentes circunferenciais

são

calculados usando-se os componentes meridionais

ainda não determinados

. Logo,

um procedimento iterativo teve que ser adotado. Primeiramente assume-se uma

distribuição uniforme para

(

−

)]). Com isso, alguns valores de

podem ser calculados usando-se as relações de grade (geometria e correlação de desvio)

estações radiais (ver item 4.4). Uma distribuição parabólica na forma

é então ajustada a esses

valores usando-se mínimos quadrados (Silva e

Miyazima, 2000). A escolha de uma distribuição parabólica é realmente adequada para

turbinas hidráulicas axiais, pois reproduz muito bem os padrões típicos de giro (Peng

., 2002a), pode recuperar o vórtice potencial (caso particular importante) e tem

fornecido boa precisão no problema em questão. Em seguida, essa distribuição

parabólica é usada no cálculo analítico de

(

) em (4.37), o que possibilita a

determinação do perfil de

(4.38) após a resolução da Equação (4.41) para

. Essa

nova distribuição de

é agora usada para recalcular os valores de

nas

grades e

as iterações continuam até que o campo de velocidade na saída do distribuidor convirja,

Figura 17.

Figura 17 – Esquema iterativo para o cálculo dos perfis de velocidade após o estator.

4.3 EQUILÍBRIO RADIAL APÓS O ROTOR

O campo de velocidade na entrada do rotor é admitido igual ao da saída do

distribuidor, aproximação razoável para a turbina considerada (Figura 11). Os perfis de

velocidade na saída do rotor são obtidos seguindo-se um procedimento análogo ao do

item 4.2 anterior. Agora, a equação da energia é escrita como (Fox e McDonald, 1998):

páSS

YYpp

++=

Pr54

(4.42)

onde

é a perda de energia mecânica por unidade de massa no rotor e

pá

é o trabalho

específico das pás, calculado usando-se a equação de Euler das turbomáquinas – a

equação fundamental das turbinas (Quantz, 1976):

)(

54 uupá

ccuY

−

(4.43)

Como fôra negligenciado o efeito das perdas hidráulicas sobre o equilíbrio radial após o

estator (Equação 4.35), resulta que a pressão de estagnação é constante na entrada do

rotor (veja a Equação 4.26):

=−==

PdSS

(4.44)

Distribuição inicial de c

Geometria da grade,

correlação de desvio;

ajuste parabólico para r c

Integral Id

c de acordo com a

condição de equilíbrio radial

e continuidade geral

Os perfis de velocidade

convergiram?

Campo de velocidade

na saída do estator

Sim

Não

Da definição de pressão de estagnação, tem-se para a saída do rotor:

S 5

ρρ

++=

(4.45)

Substituindo a equação de equilíbrio radial, dada por:

(4.46)

e as Equação (4.44) e (4.45) na derivada da equação da energia (4.42), obtém-se:

pá

ρρρρρ

++++=

(4.47)

Da Equação (4.43), deriva-se:

ucd

pá

)()(

−=

(4.48)

O primeiro termo do segundo membro de (4.48) já está disponível, pois a distribuição

é igual a de

, que fora calculada no item 4.2 anterior usando-se um ajuste

parabólico (

). Substituindo então (4.48) em

(4.47), resulta:

ucd

rcc

)()(

+++−+=

(4.49)

equação essa que expressa a condição de equilíbrio radial após o rotor. Novamente,

distribuições de

que satisfaçam a essa equação estarão em conformidade com

os princípios de conservação de energia e quantidade de movimento, constituindo um

campo de velocidade realista para o escoamento na saída do rotor da turbina. Pode-se

ainda rearranjar (4.49) na forma:

ucd

Ycd

rcd

umu

)()2/()(

)(0

4Pr

+−=

(4.50)

ucd

rcd

Ycd

uuum

)()()()2/(

455Pr

−

−=

(4.51)

que é mais conveniente para se integrar e expressar

em termos de

. Novamente,

pelas mesmas justificativas dadas no item 4.2, o efeito das perdas hidráulicas sobre o

equilíbrio radial após o rotor pode ser negligenciado, isto é, admite-se

/ dr

= 0 na

Equação (4.51). Mais uma vez, não se está afirmando serem nulas as perdas no rotor,

mas apenas que seu efeito sobre o equilíbrio radial é desprezível (Peng

et al

., 2002a).

Dessa forma, tem-se após a integração do cubo (

) a um raio

genérico:

()

ruru

hmm

ucucdr

rcd

crc ||2

)()(

2)(

−−

−

−=−

∫

(4.52)

Definindo:

()

ruru

ucucdr

rcd

rI ||2

)()(

2)(

−−

−

−=

∫

(4.53)

escreve-se explicitamente a distribuição de velocidade meridional após o rotor

satisfazendo ao equilíbrio radial ao longo de todo o raio:

)()(

rIcrc

rhmm

(4.54)

A continuidade geral é novamente imposta para o cálculo da velocidade meridional na

estação do cubo (

), levando à Equação não-linear (4.57), análoga à obtida para o

distribuidor (4.41):

∫

rdrrcQ )(2

(4.55)

)(

rdrrc

∫

(4.56)

)(

rdrrIc

rhm

=⋅+

∫

(4.57)

A Equação (4.57) também é resolvida usando-se o método de bisseção da função

fzero

do MatLab

e a integração em (4.57) é realizada numericamente segundo a

regra de Simpson.

Analogamente ao cálculo da integral

(

) em (4.37), o cálculo da integral

(

)

em (4.53) apresenta a mesma dificuldade inicial: o conhecimento prévio da distribuição

do momento cinemático

(

). Porém, como ocorre na saída do distribuidor, os

componentes circunferenciais (

) são calculados usando-se os componentes

meridionais (

)

ainda não determinados

. Adotou-se, então, o mesmo esquema

iterativo para resolver o problema. Primeiramente, assume-se uma distribuição inicial

para

(por exemplo,

). Com isso, alguns valores de

podem ser

calculados usando-se as relações de grade (geometria, triângulos de velocidade e

correlação de desvio) em

estações radiais (ver 4.4). Uma distribuição cúbica na forma

é então ajustada a esses

valores usando-se mínimos

quadrados. A escolha de uma distribuição cúbica em vez de uma parabólica mostrou-se

necessária para reproduzir as inflexões típicas no perfil de

. Em seguida, essa

distribuição cúbica é usada no cálculo analítico de

(

) em (4.53), onde a regressão

parabólica para

) também é utilizada. Isso possibilita a determinação do

perfil de

(4.54) após a resolução da Equação (4.57) para

. Essa nova distribuição

é agora usada para recalcular os valores de

nas

grades e as iterações

continuam até que o campo de velocidade na saída do rotor convirja, Figura 18.

Figura 18 – Esquema iterativo para o cálculo dos perfis de velocidade após o rotor.

Para o rotor, porém, um esquema iterativo puramente de substituição apresentou

sérios problemas de convergência. Isso se deve ao caráter altamente não-linear das

equações e ao fato de que o perfil de velocidade na entrada do rotor é não-uniforme,

diferentemente do que ocorre na entrada do estator, havendo ainda absorção da energia

do fluido pelas pás do rotor. A saída escolhida para contornar essa dificuldade, que

surgiu devido ao modo como se resolveu a condição de equilíbrio radial, foi empregar

um esquema de sub-relaxação. Embora isso retarde a convergência do algoritmo, evitou

as instabilidades que a impediam. Cada vez que uma nova distribuição de

Distribuição inicial de c

Geometria da grade,

triângulos de velocidade,

correlação de desvio;

esquema de sub-relaxação,

ajuste cúbico para

r c

Integral Ir

c de acordo com a

condição de equilíbrio radial

e continuidade geral

Os perfis de velocidade

convergiram?

Campo de velocidade

na saída do rotor

Sim

Não

calculada (Equação 4.54), levando a uma nova distribuição de

nas grades (a partir

dos triângulos de velocidade;

grade

), os novos valores atribuídos de fato para

são:

atual

grade

novo

ccc

555

)1(

λλ

−+=

(4.58)

onde

é o fator de sub-relaxação. Para iniciar esse esquema, a primeira distribuição de

é igualada a zero. Como faixa recomendada para os valores do fator de sub-

relaxação, tem-se 0,05

≤

0,20, sendo que 0,10 foi o valor usado nos resultados a

serem apresentados.

4.4 CORRELAÇÃO DE DESVIO; TRIÂNGULOS DE

VELOCIDADE

Os desvios do escoamento em relação aos ângulos geométricos de saída das

grades, tanto no distribuidor como no rotor, são avaliados usando-se a correlação de

desvio de Carter e Hughes (Horlock, 1973). Por exemplo, em uma dada estação radial

do distribuidor, o ângulo de desvio na saída,

(Figura 19), é dado por:

aletaf

φααδ

=−=

(4.59)

Figura 19 – Desvio angular do escoamento na saída de uma grade do distribuidor.

aleta

onde

é o ângulo do escoamento,

aleta

é o ângulo geométrico (ângulo construtivo)

da palheta diretriz,

aleta

−

aleta

)

é o ângulo de arqueamento do perfil (

camber

angle

é a corda do perfil,

é o passo da estação radial e

é um fator empírico. Em

Horlock (1973),

é fornecido graficamente em função do ângulo de montagem

(em

graus;

stagger angle

) e do tipo de linha média do perfil (circular ou parabólica). Esse

gráfico foi aproximado pela função linear

(

) = 0,21 – 0,04(90 –

)/60, para as linhas

de arqueamento ARC adotadas neste trabalho (capítulo 3).

Uma vez determinado

aleta

) e tendo-se a distribuição de

(item

4.2), calcula-se (veja a Figura 14):

sen

(4.60)

tan

(4.61)

Com isso, tem-se também a velocidade absoluta na entrada do rotor, que foi admitida

igual à velocidade na saída do distribuidor (item 4.3), ou seja:

⎪

⎩

⎪

⎨

⎧

⇔=

αα

(4.62)

Daí, calcula-se imediatamente (veja a Figura 15):

⎟

⎠

⎞

⎜

⎝

⎛

−

arctan

(4.63)

sen

(4.64)

O ângulo do escoamento relativo na saída do rotor é determinado também por

meio da correlação de desvio de Carter e Hughes (Figura 20):

555

pápáf

φβδββ

+=+=

(4.65)

onde agora

pá

−

pá

(

) = 0,21 – 0,04(90 –

)/60, sendo

o ângulo de

montagem da pá (em graus) e

pá

, respectivamente, seus ângulos geométricos de

entrada e de saída. Tendo-se também a distribuição de

(item 4.3), calcula-se (veja a

Figura 15):

sen

(4.66)

tan

−=

(4.67)

55 mu

ccc +=

(4.68)

⎟

⎠

⎞

⎜

⎝

⎛

−=

arctan2/

πα

(4.69)

finalizando os elementos dos triângulos de velocidade necessários aos cálculos de

deflexões, incidências, ângulos de ataque, trabalhos específicos, perdas, etc.

Figura 20 – Desvio angular do escoamento na saída de uma grade do rotor.

pá

Deve-se comentar ainda que a correlação de Carter e Hughes foi desenvolvida

sobre dados experimentais de turbinas a gás subsônicas (de baixa velocidade) e, como

não há dependência explícita dos parâmetros dessa correlação com o número de Mach

ou qualquer outra grandeza relacionada à compressibilidade do escoamento, admitiu-se

ser aplicável também a escoamentos incompressíveis (turbinas hidráulicas). Se algum

parâmetro dependesse, por exemplo, do número de Mach, então o número de Mach

seria igualado a zero numa tentativa de ainda assim aplicar essa correlação. De fato,

modelos de desvios e de perdas específicos para turbinas hidráulicas são bastante

escassos na literatura em comparação com as numerosas publicações referentes a

turbinas a gás.

Saliente-se também que o modelo de desvio de Carter e Hughes só é válido em

situações afastadas do

stall

da grade em questão, ou seja, em situações relativamente

próximas de incidência zero e com deflexões moderadas. Note também que, nesse

modelo, o desvio depende somente da geometria da grade (passo, corda, ângulo de

montagem, ângulo de arqueamento, tipo de linha de esqueleto), o que torna mais direta

sua aplicação.

4.5 CORRELAÇÕES DE PERDAS; RENDIMENTO DA

TURBINA

As perdas hidráulicas ou de escoamento através da turbina são avaliadas usando-

se correlações empíricas da literatura. Há várias formas para se classificar e tratar as

perdas em turbinas hidráulicas (Ueda, 1982). Procurou-se adotar um modelo que

contemple as parcelas mais significativas em turbinas axiais do tipo tubular (Figura 11).

Para o distribuidor, admite-se apenas perdas de perfil, causadas pelo atrito viscoso com

as aletas. Já para o rotor, além das perdas de perfil, consideram-se as perdas por choque

na entrada, ocasionadas por ângulos desfavoráveis ao escoamento incidente sobre as

pás. Tem-se também as perdas no tubo de sucção, associadas à recuperação imperfeita

de pressão estática nesse difusor e que, em turbinas de baixa queda, representam

normalmente a maior parcela da perda total (Drtina e Sallaberger, 1999; Moura e Brasil

Júnior, 2005). As perdas são calculadas após a determinação dos perfis de velocidade

nas saídas do distribuidor e do rotor, pois são correlacionadas às velocidades nessas

seções.

Ao se usarem modelos de perdas em um projeto, o fundamental é que esses

modelos indiquem as corretas tendências de variação das perdas causadas por alterações

na geometria da turbina. O conjunto de correlações escolhido deve orientar a busca em

direção a uma geometria de projeto de boas características, ou seja, com perdas

mínimas. O modo como cada parâmetro de projeto afeta as perdas deve ser evidenciado

pelo modelo adotado. Por exemplo, sabe-se que as perdas no tubo de sucção devem-se,

em grande parte, ao giro do escoamento em sua entrada. Assim, um modelo de perdas

para o tubo de sucção deve conduzir, durante o processo de otimização, a uma

geometria do rotor que produza pouco giro na saída (pequenos valores para

próximo de 90°). Da mesma forma, o modelo de perda por choque na entrada do rotor

deve levar a uma geometria para as pás com ângulos de montagem e arqueamentos tais

que se produzam simultaneamente pequenos valores de incidência na entrada e

deflexões resultantes ainda assim razoáveis (para que haja bons valores de trabalho

específico). Citando Denton (1993): “...o uso de correlações não deve ser um substituto

de se tentar compreender as origens das perdas... uma boa compreensão física para elas

pode ser mais valiosa do que uma predição quantitativa.” Também, “é importante ao

usar as correlações, reconhecer suas limitações e não desenvolver um falso senso de

segurança se elas fornecem respostas corretas...” O mais importante é “...ser aptos a

identificar características boas e ruins do escoamento e alterar nossos projetos de

acordo, mesmo que não possamos quantificar as melhorias antes de testá-las.”

Finalmente, “...não devemos temer admitir nossa falta de compreensão dos mecanismos

complexos de geração de entropia. O progresso ocorre mais facilmente quando estamos

cientes de nossas limitações e continuamente nos esforçamos para reduzi-las.”

A seguir, descrevem-se as correlações empregadas no cálculo das perdas

específicas (J/kg) nos componentes da turbina. Depois, explica-se a forma como essas

perdas são integradas nos fluxos de massa para se obter a potência perdida, em watts, e

o modo como se calcula a potência absorvida pelo rotor. Feito isso, determina-se a

altura de energia da turbina (trabalho específico da máquina) e seu rendimento total.

4.5.1 – Perda no distribuidor

Apenas perdas por atrito (perdas de perfil) são consideradas para o distribuidor.

Ignoram-se, por exemplo, perdas por choque na entrada, admitindo-se implicitamente

que o fluxo incida de maneira ideal sobre as palhetas diretrizes, hipótese razoável para

turbinas hidráulicas axiais (geralmente com apenas um único estágio).

A perda por unidade de massa em uma estação radial,

, é dada por:

dPd

(4.70)

onde

é a velocidade absoluta na saída da grade e

é um coeficiente empírico. No

caso,

é calculado usando-se a correlação de Soderberg (Horlock, 1973), descrita mais

adiante.

4.5.2 – Perda no rotor

Para o rotor, além das perdas por atrito, é essencial considerar também as perdas

por choque na entrada, associadas a uma incidência desfavorável do escoamento sobre

as pás. De fato, se a perda por choque não for levada em consideração, geometrias com

altas incidências no rotor e separação precoce do escoamento poderiam não ser

descartadas durante o processo de otimização, comprometendo a otimalidade prática do

projeto encontrado.

Para a perda de perfil (

) o cálculo é análogo ao do distribuidor. Em uma dada

estação, tem-se:

(4.71)

onde

é a velocidade relativa na saída da grade considerada e

é calculado usando-

se a correlação de Soderberg.

Se o número de pás fosse infinito, a condição sem choque seria aquela em que a

velocidade relativa na entrada do rotor fosse tangente à pá, ou seja, com

pá

(incidência zero). Com um número finito de pás, a condição sem choque passa a ocorrer

com incidências levemente diferentes de zero e depende essencialmente do número de

pás e da forma dos perfis. Para o nível de otimização proposto neste trabalho, será

suficiente aproximar a condição sem choque daquela com incidência nula, embora não

a seja exatamente. O importante aqui é descartar geometrias com incidências realmente

muito desfavoráveis, ainda que a incidência ótima não possa ser determinada com

precisão. Posto isso, o cálculo da perda por choque,

Pch

, é realizado usando-se o

esquema proposto por Pfleiderer e Petermann (1979), mas com o componente de

choque,

, sendo considerado igual ao componente de incidência:

Pch

(4.72)

0,5

≤

0,7 e:

pá

incch

−

⎟

⎠

⎞

⎜

⎝

⎛

+==

tantan

αβ

(4.73)

4.5.3 – Perda no tubo de sucção

O escoamento em tubos de sucção é complexo, com importantes efeitos

rotacionais tridimensionais e pulsações geradas por instabilidades das estruturas

vorticais. A avaliação precisa das perdas no tubo de sucção normalmente exige o uso de

métodos CFD com modelos de turbulência bem aferidos (Moura e Brasil Júnior, 2005).

Para os propósitos de projeto do trabalho presente, no entanto, será satisfatório o

procedimento proposto por Raabe (1985), no qual as perdas no tubo de sucção,

Pts

, são

associadas de uma forma bastante simples aos componentes circunferencial e

meridional da velocidade absoluta na saída do rotor, numa tentativa de se considerar as

contribuições do giro e da difusão do escoamento (Ueda, 1982). Numa dada estação

radial, tem-se:

5 u

DmPts

XY +=

(4.74)

sendo

, respectivamente, os componentes meridional e circunferencial da

velocidade absoluta na saída da grade do rotor.

é o coeficiente de perda por difusão

e depende basicamente da geometria do tubo de sucção, podendo ser adotado em torno

de 0,09 para tubos de sucção retilíneos (típicos nas turbinas dos tipos bulbo e tubular)

ou 0,12 para tubos de sucção em cotovelo (típicos nas turbinas axiais convencionais de

eixo vertical).

é o coeficiente de perdas por giro, recomendado entre 0,2 e 0,4.

4.5.4 – Perda mecânica

A perda mecânica ou externa, causada pelos atritos nos mancais, vedações e

acoplamentos da turbina, é considerada simplesmente admitindo-se de antemão um

valor para o rendimento mecânico,

mec

0,99 a 0,96=

mec

(4.75)

4.5.5 – Correlação de Soderberg

A correlação de Soderberg (Horlock, 1973) é usada neste trabalho para estimar

as perdas por atrito, tanto no distribuidor como no rotor. Embora tenha sido

desenvolvida para aplicações em turbinas a gás (assim como a correlação de desvio de

Carter e Hughes, item 4.4), ela se baseia em parâmetros de carregamento de grades

genéricas (tipo injetoras), com razão de aspecto finito e sem envolver explicitamente

qualquer parâmetro de compressibilidade, podendo servir a aplicações também em

turbinas hidráulicas. Cite-se novamente a relativa escassez de publicações referentes a

modelos de perdas específicos para turbinas hidráulicas, o que levou ao uso de

correlações oriundas de testes com turbinas a gás.

Os efeitos do número de Reynolds e da razão de aspecto das grades (

; razão

entre as alturas radial e axial de uma pá) são levados em conta na correlação de

Soderberg. O número de Reynolds,

, é calculado com base num diâmetro hidráulico

e numa velocidade

referente à seção de saída das grades:

R =

(4.76)

sendo

a viscosidade dinâmica do fluido e

, dependendo de se tratar do

distribuidor ou do rotor, e:

cos

cos2

cos

cos2

(4.77)

sendo

o passo da grade em questão.

A correlação completa de Soderberg parte inicialmente de uma sub-correlação

para o coeficiente de perdas em grades padronizadas com razão de aspecto 3:1,

operando com

= 10

. Para essas grades, o coeficiente de perdas,

, depende

fortemente da deflexão angular do escoamento através da grade,

−

, no

distribuidor;

–

, no rotor), e depende fracamente da espessura dos perfis. Essa

sub-correlação está apresentada graficamente em Horlock (1973). Uma expressão

aproximada foi obtida para representar esse gráfico:

ξξ

01053,0

(4.78)

onde

está em graus. O coeficiente de perdas para deflexão nula,

, pode ser ajustado

em função do grau de aperfeiçoamento tecnológico da turbina (quanto menor o

coeficiente, maior o grau de aperfeiçoamento). Admite-se 0,04

≤

0,06.

A correlação completa de Soderberg avalia a perda em uma grade,

, com

qualquer relação de aspecto

B/b

, operando com qualquer número de Reynolds:

⎥

⎦

⎤

⎢

⎣

⎡

−

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

= 1075,0975,0)1(

4/1

ξξ

(4.79)

A Tabela 1 a seguir reúne todas as correlações empíricas de perdas usadas neste

trabalho.

Tabela 1 – Conjunto de modelos de perdas para a turbina

Mecanismo de perda Correlação Referência

Perda de perfil no estator

(atrito nas aletas)

dPd

Horlock (1973)

Perda por choque no rotor

(incidência nas pás)

Pch

onde

= 0,5 a 0,7 e

pá

incch

−

⎟

⎠

⎞

⎜

⎝

⎛

+==

tantan

αβ

Pfleiderer e Petermann

(1979)

Perda de perfil no rotor

(atrito nas pás)

Horlock (1973)

Perda no tubo de sucção

(difusão e giro)

5 u

DmPts

XY +=

com

= 0,09 ou 0,12 e

= 0,20 a 0,40

Raabe (1985)

Perda mecânica

(perda externa)

0,99 a 0,96

mec

Adotado

Coeficientes da correlação de Soderberg:

⎥

⎦

⎤

⎢

⎣

⎡

−

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

= 1075,0975,0)1(

4/1

ξξ

01053,0

= 0,04 a 0,06

−

(deflexão)

= razão de aspecto da grade (radial/axial)

= viscosidade

cos

cos2

cos

cos2

Horlock (1973)

4.5.6 – Integração das perdas e da potência absorvida pelas pás

As perdas hidráulicas por unidade de massa descritas no item anterior são

calculadas em

estações radiais (

grades) a fim de serem integradas nos fluxos de

massa para se obter a potência perdida total na turbina,

. Como as perdas no rotor e

no tubo de sucção são calculadas no mesmo fluxo de massa (seção de saída do rotor),

elas são somadas. Usando-se então mínimos quadrados, ajusta-se um polinômio do

terceiro grau aos

(

≥

4) valores de perda específica no distribuidor e outro polinômio

do terceiro grau aos

(

≥

4) valores de perda específica no (rotor + tubo de sucção). A

escolha por simples funções cúbicas mostrou-se bastante satisfatória nesses ajustes.

Outra forma de se computar a potência perdida é calcular as perdas específicas em um

número razoavelmente grande de estações e usar esses valores diretamente num

esquema de integração numérica. Porém, fazer uso de regressões para bem representar a

distribuição radial das perdas possibilita sua avaliação em apenas poucas estações, o

que é importante quando esses cálculos demandam grande esforço computacional

(embora este não seja o caso).

Realizadas as regressões cúbicas para

e (

Pch

Pts

), essas funções são

usadas nas integrações das perdas:

∫∫

mPdPdPd

rdrcYmdYP

2 Seção

πρ

(4.80)

∫∫

++=++=

mPtsPchPtsPchtsrP

rdrcYYYmdYYYP

5Pr

5 Seção

Pr)(

)(2)(

πρ

(4.81)

As velocidades meridionais são calculadas usando-se as Equações (4.38) e (4.54), que

fornecem as distribuições de

, respectivamente, de acordo com a condição de

equilíbrio radial e satisfazendo à continuidade geral (veja os itens 4.2 e 4.3). As

integrações em (4.80) e (4.81) são feitas numericamente via regra de Simpson.

A potência perdida total,

, é então dada por:

)( tsrPPdP

PPP

(4.82)

A potência absorvida pelo rotor,

pá

, resulta da integração dos trabalhos

específicos

pá

(Equação 4.43). Tem-se:

∫∫∫

−==

5 Seção

4 Seção

Rotor

mducmducmdYP

uupápá

&&&

(4.83)

∫∫

−=

mupá

rdrcucrdrcucP

5544

πρπρ

(4.84)

Os valores de

) e

são calculados usando-se as regressões polinomiais

determinadas para

) e

), descritas

respectivamente nos itens 4.2 e 4.3. Os componentes meridionais são novamente dados

pelas Equações (4.38) e (4.54).

Uma maneira aproximada de calcular

pá

é:

∫

−=

muupá

rdrcccuP )(2

πρ

(4.85)

com

2/)(

54 mmm

ccc +=

. Esse procedimento equivale à Equação (4.83) quando as

linhas de corrente são paralelas no plano meridional, o que levaria a

(continuidade nos tubos de corrente). Em turbinas hidráulicas axiais, as linhas de

corrente não são exatamente paralelas no plano meridional, mas apresentam curvaturas

realmente suaves (Drtina e Sallaberger, 1999; Osterwalder, 1960), de modo que as

distribuições de

não são muito diferentes e a Equação (4.85) pode ser aplicada

com precisão razoável.

4.5.7 – Rendimento e altura de energia da turbina

Uma vez determinadas a potência absorvida pelo rotor (

pá

) e as perdas

hidráulicas totais (

), calcula-se a potência útil ou hidráulica,

, por:

Ppáh

PPP +=

(4.86)

A potência útil ou hidráulica é a taxa de extração da energia mecânica do fluido em

escoamento através da turbina. Representa a energia de “entrada” da máquina. A fração

dessa energia que é efetivamente absorvida pelas pás do rotor é o chamado rendimento

hidráulico,

, um importante parâmetro de eficiência da turbina:

pá

(4.87)

De fato, o projetista busca uma configuração geométrica para as grades do distribuidor

e do rotor, assim como uma forma para o tubo de sucção, de tal modo que se obtenha o

máximo rendimento hidráulico, isto é, a máxima absorção da energia mecânica

disponível no escoamento pelas pás do rotor.

O rendimento total da turbina,

, é fração da potência hidráulica transformada

em potência de eixo,

, a potência de “saída” da máquina. Pode ser calculado por:

mech

ηηη

(4.88)

sendo

mec

a eficiência mecânica ou externa da turbina.

O rendimento total constitui a função objetivo do presente problema de

otimização, formulado no capítulo 2 (na verdade, o valor objetivo é –

, pois o problema

foi posto na forma de minimização). Como

mec

é previamente fixado (veja o item

4.5.4), a maximização de

equivale à maximização de

. De fato, esse é o objetivo

maior do projeto hidrodinâmico de uma turbina.

Deve-se comentar que as perdas por fuga e por atrito lateral não são

consideradas aqui, embora sejam típicas no estudo de turbomáquinas e desempenhem

papel importante em seu desempenho (Bran e Souza, 1979; Denton, 1993). Isso porque,

numa fase inicial de projeto, o contexto deste trabalho, a preocupação maior de fato é

com relação às perdas de escoamento, que determinam o rendimento hidráulico. A

minimização das perdas por fuga e por atrito lateral freqüentemente está ligada mais a

questões de fabricação, montagem e manutenção do que ao projeto hidrodinâmico

propriamente dito, fugindo do escopo deste trabalho.

O trabalho específico da turbina,

, que constitui a energia mecânica por unidade

de massa que o fluido disponibiliza à máquina de fluxo, é definido por:

(4.89)

Quando expresso como energia por unidade de peso (J/N ou m), denomina-se altura de

energia, ou altura de queda disponível, ou simplesmente altura,

(4.90)

sendo

a aceleração da gravidade local. A altura

está sujeita a restrições de

desigualdade no problema de otimização (veja as Equações 2.1 a 2.4).

A seqüência geral dos cálculos desenvolvidos neste capítulo para se obter o

rendimento e a altura da turbina, a partir das variáveis de projeto e demais informações,

é resumida no fluxograma da Figura 21.

Figura 21 – Esquema geral da seqüência de cálculos para a análise do escoamento

(para os perfis de velocidade, há ainda os esquemas iterativos das

Figuras 17 e 18).

Variáveis de projeto (xj)

Dados: variáveis de projeto, outras variáveis,

propriedades, fatores empíricos,

Q, n

Parametrização da geometria das pás

Perfis de velocidade na saída do distribuidor

Correlação

de desvio

Perfis de velocidade na saída do rotor

Integração das perdas e

dos trabalhos específicos

Correlações

de perdas

Rendimento (

η)

Altura (

Cálculo do

escoamento

Capítulo 5

MÉTODOS DE OTIMIZAÇÃO

Os algoritmos de otimização são desenvolvidos para substituir uma parte onerosa

do trabalho do projetista. Constituem a forma automática de se alterar as variáveis de

projeto, analisar as conseqüências dessas alterações e tornar a modificá-las no intuito de

se cumprir o objetivo do problema, no caso, maximizar o rendimento da turbina. Neste

capítulo, são apresentados os dois métodos de otimização usados neste trabalho, a

saber, um Algoritmo de Programação Quadrática Seqüencial (SQP) e um Algoritmo de

Busca Aleatória Controlada (CRSA).

5.1 INTRODUÇÃO

Métodos de otimização local (baseados em gradiente) estão bem desenvolvidos

atualmente. Há diversos resultados matemáticos que garantem a existência de ótimos

locais e a convergência dos métodos de busca local (Nash e Sofer, 1996). Diversos

pacotes comerciais de

softwares

matemáticos, como o Excel



, MatLab



Mathematica



, apresentam

toolboxes

de otimização, seja para programação linear ou

não-linear, com ou sem restrições, mas quase todos destinados apenas a buscas locais.

Em muitas aplicações, o uso de métodos locais pode ser satisfatório. Por

exemplo, quando é bem conhecido o comportamento da função objetivo e das restrições

do problema, pode-se iniciar o método de otimização local em um ponto que já se sabe

estar próximo da solução. Com um número não tão grande de tentativas, às vezes é

realmente possível investigar o ótimo global procurando-se por ótimos locais. Em

outras situações, a preocupação maior pode não ser o ótimo global propriamente dito,

mas uma solução mais adequada que combine um valor razoável de função objetivo e

uma menor variabilidade em sua vizinhança. Tal consideração é mais relevante quando

as variáveis do problema representam grandezas físicas ou itens de produção,

naturalmente associados a valores médios e distribuições de probabilidade (devido por

exemplo a não homogeneidade dos materiais, tolerâncias de fabricação e montagem,

oscilações na oferta e demanda de produtos, etc.). Esses aspectos podem ser

considerados no contexto da

otimização robusta

Por outro lado, em alguns problemas é desejável e útil encontrar realmente um

ótimo global. Esse é o caso, por exemplo, do projeto de turbinas hidráulicas, onde

pequenos aumentos de rendimento podem representar grandes aumentos na energia

produzida. De fato, um dos propósitos deste capítulo é apresentar as vantagens em usar

um algoritmo de busca global em comparação com as buscas locais no problema de

otimização do projeto de turbinas hidráulicas axiais.

O algoritmo de busca local escolhido corresponde à função

fmincon

toolbox

de otimização do MatLab



(versão 5.3). Consiste em uma programação quadrática

seqüencial com restrições, usando-se busca em linha, calculando-se as derivadas

direcionais por diferenças finitas e aproximando-se a matriz Hessiana pela fórmula de

BFGS (Nash e Sofer, 1996). As duas principais desvantagens dessa escolha pelo

método de otimização são a busca apenas por mínimos locais e a necessidade de se

fornecer um ponto de partida. Essas duas características, em conjunto, tornam o êxito

do processo de otimização bastante dependente da estimativa inicial para as variáveis

de projeto. O projetista, assim, deve conhecer de antemão uma configuração não muito

afastada da ótima. Também, a solução obtida, presumindo-se a convergência do

algoritmo, pode não ser o ótimo global do problema. Além disso, investigações

preliminares no problema da turbina axial mostraram que mesmo pontos de partida

ligeiramente diferentes podem levar a soluções distintas, com distintos valores de

rendimento. Em outras palavras, o problema da turbina hidráulica axial apresenta

diversos mínimos locais (a função objetivo, nesse caso, é multimodal, “rugosa”).

Para superar essas limitações, propõe-se também o uso de um método de busca

global. Entre as possibilidades disponíveis, escolheu-se um Algoritmo de Busca

Aleatória Controlada (

Controlled Random Search Algorithm

– CRSA), a saber, o CRSI

de Ali

et al

. (1997) com algumas modificações. Apesar de ser uma opção menos

conhecida na comunidade de engenharia, o CRSA foi escolhido em virtude de sua

facilidade de implementação, relativa rapidez e bons resultados reportados na literatura

(Ali e Törn, 2004; Manzanares Filho

et al

., 2005). CRSA’s são algoritmos

populacionais, como os algoritmos Genético e de Evolução Diferencial. Neles, parte-se

de uma população inicial de pontos sobre a região viável do problema e então procede-

se a substituições iterativas dos piores pontos por pontos tentativa melhores, na

esperança de que a população contraia-se em torno do ótimo global.

Como ocorre também nos outros algoritmos evolucionários, um CRSA não

possui garantias de convergência ou de otimalidade (Ali

et al

., 1997). De fato, um dos

poucos resultados para ótimos globais é o seguinte teorema (Nash e Sofer, 1996):

Soluções globais de programas convexos. Seja

um mínimo local de um problema de

programação convexa. Então

é também um mínimo global. Se a função objetivo é

estritamente convexa, então

é o único mínimo global

. Esse teorema é importante em

programação linear (caso particular de programação convexa); no problema de

otimização da turbina hidráulica axial, porém, a função objetivo não é provavelmente

nem côncava nem convexa e esse teorema não poderia ser aplicado.

Apesar da natureza heurística do CRSA, diversos estudos têm demonstrado a

eficácia do algoritmo em buscas globais e sua capacidade de aplicação em problemas

complexos de engenharia (Ali

et al

., 1997, Manzanares Filho

et al

., 2005).

A principal vantagem do CRSA é que o método tem partida automática. Em vez

de uma estimativa inicial (que o projetista deveria fornecer com um certo cuidado), o

algoritmo usa uma população inicial aleatoriamente distribuída sobre a região de busca

do problema. Do ponto de vista do engenheiro, escolher faixas para as variáveis de

projeto é mais fácil do que determinar uma boa estimativa inicial. Além disso, há a

esperança de se encontrar uma solução global, mesmo com regiões de busca levemente

diferentes e diferentes populações iniciais. Um número de chamadas da função objetivo

relativamente pequeno também é uma expectativa quando se usa CRSA (Ali e Törn,

2004). Isso é particularmente importante quando a função objetivo demanda um alto

esforço computacional – algo que freqüentemente ocorre na otimização de

turbomáquinas.

5.2 PROGRAMAÇÃO QUADRÁTICA SEQÜENCIAL

Sequential Quadratic Programming

(SQP) é um dos métodos mais difundidos

para a busca por ótimos locais em problemas de minimização não-linear com restrições,

como é o caso do problema da turbina hidráulica axial estabelecido no capítulo 2. O

SQP foi sugerido inicialmente por Wilson (1963) e desde então diversos autores têm

dado contribuições. Enquadra-se na categoria de Métodos de Pontos Viáveis (

Feasible-

Point Methods

), entre os quais inclui-se também os Métodos de Gradiente Reduzido

(

Reduced-Gradient Methods

). A idéia principal do SQP é obter uma direção de busca

resolvendo-se um programa quadrático, i.e., um problema com função objetivo

quadrática e restrições lineares, constituindo uma generalização dos métodos quase-

Newton para minimização sem restrições (Nash e Sofer, 1996).

As linhas gerais do método descrevem-se a seguir (adaptado de Sousa Júnior

., 2005). O problema é formulado analogamente à Expressão (2.1), incluindo-se agora

restrições de igualdade:

minimizar

(

)

sujeito a

(

) = 0 (5.1)

(

)





é o vetor

-dimensional das variáveis,

(

) é a função objetivo,

(

) e

(

) são os

vetores de funções de restrição e





é a região de busca do problema, definida por

restrições laterais.

O SQP baseia-se em obter subproblemas usando-se uma aproximação quadrática

do lagrangiano da função objetivo e linearizando-se as restrições, ou seja:

minimizar

pxfpBp

)(



sujeito a

0)()(



xgpxg

(5.2)

0)()(



xgpxg

kukl

xxpxx





onde

é uma aproximação positiva semi-definida da matriz Hessiana,

é a iteração

atual e

são vetores contendo os limites inferior e superior das variáveis.

Seja

a solução do subproblema (5.2). Uma busca uni-direcional (

line search

) é

usada na procura de um novo ponto

)1;0(







kkkk

pxx

(5.3)

tal que uma “função mérito” (

merit function

) tenha um valor mínimo nesse novo ponto.

A função de Lagrange aumentada é comumente usada como

merit function

. Enquanto a

otimalidade não for obtida, as aproximações

para a matriz Hessiana são atualizadas

segundo a fórmula de BFGS (Nash e Sofer, 1996).

O MatLab



dispõe da função

fmincon

em seu

toolbox

de otimização.

fmincon

implementa um método SQP na forma descrita acima e foi considerada a função mais

conveniente para as buscas locais neste trabalho.

fmincon

busca um mínimo local de

uma função escalar de várias variáveis, sujeita a restrições, partindo-se de uma

estimativa inicial

. Constitui o bloco central do processo de otimização, alterando as

variáveis de projeto

, enviando-as ao programa de análise do escoamento (o

solver

Figura 21), avaliando o rendimento alcançado e tornando a modificar a geometria em

busca de uma configuração ótima, respeitando-se as restrições laterais para as variáveis

de projeto e as restrições de desigualdade para a altura de energia, Figura 22.

Figura 22 – Fluxograma do projeto otimizado usando-se SQP-

fmincon

Estimativa inicial x ,

faixas para as variáveis de projeto (restrições

laterais), opções para a função

fmincon

Otimizador (fmincon)

Ciclo interno de convergência

(SQP) em direção a um ótimo da

função objetivo

Função objetivo

solver

Projeto otimizado,

relatório diagnóstico do processo de otimização

Restrições não-lineares

Variáveis de

projeto,

Rendimento negativo



Vetor de funções

de restrição,

Convergência,

solução

Neste trabalho, optou-se pelo cálculo numérico dos gradientes da função

objetivo, usando-se um esquema de diferenças finitas. Isso constitui o que se chama de

otimização de “média escala” no

fmincon

. Maiores detalhes sobre o algoritmo são

encontrados no manual do

Optimization Toolbox

do MATLAB



Saliente-se novamente que as duas principais desvantagens dessa técnica de

otimização são justamente a busca apenas por mínimos locais e a necessidade de se

fornecer um ponto de partida. Em conjunto, essas características tornam o êxito do

processo de otimização bastante dependente da estimativa inicial

para as variáveis de

projeto. Isso motivou o estudo de um método de busca global para o problema da

turbina hidráulica axial.

5.3 ALGORITMOS DE BUSCA ALEATÓRIA

CONTROLADA (CRSA)

(Adaptado de Manzanares Filho

et al

., 2005).

Controlled Random Search Algorithms

–

CRSA – são algoritmos adequados para a otimização global de funções reais contínuas







, não necessariamente diferenciáveis, definidas no “hiper cubo” (

hyper-box

)

= {





xxx 

= 1, ...,

}, onde

representam, respectivamente,

limites inferiores e superiores para as

coordenadas de

. Um ponto

* é dito ser um

mínimo global de

(



(





. Além das restrições laterais usadas na

definição de

, outros tipos de restrição poderiam, em princípio, ser impostos por meio

de um esquema de penalização na função objetivo, como na Equação (2.4).

Os CRSA’s foram propostos como aperfeiçoamento para os métodos de busca

aleatória simples, nos quais apenas o ponto com o melhor valor de função objetivo

permanece em cada iteração (Price, 1977). Como os algoritmos Genético e de Evolução

Diferencial, um CRSA é um algoritmo populacional que parte de um conjunto ou

população inicial

pontos aleatoriamente tomados em

e então executa um

processo iterativo de contração dessa população em direção ao ótimo global por meios

puramente heurísticos (Ali

et al.

, 1997; Ali e Törn, 2004). Nos CRSA’s, o tamanho

da população é mantido ao longo de todo o processo de otimização. Diferentemente dos

outros algoritmos populacionais mencionados, um CRSA substitui um único ponto da

população (seu atual pior ponto,

, com o maior valor de função objetivo) por um ponto

melhor

em cada iteração (i.e., um ponto tentativa

tal que

(

)



(

)). Logo sua

implementação é mais direta.

5.3.1 – Algoritmo básico

O CRSA básico para minimização é descrito como se segue (adaptado de Ali

al.

, 1997, e Ali e Törn, 2004):

Gere uma população inicial

pontos aleatórios em

= {

, ...,

Calcule os valores de função objetivo nesses pontos de uma forma indexada.

Determine o pior ponto,

, e o melhor ponto,

, i.e., aqueles pontos em

com o

maior e o menor valores de função,

, respectivamente. Se um critério de

parada já estiver sendo cumprido, então pare (por exemplo, pare se





onde



é uma tolerância dada).

Gere um ponto tentativa

para substituir o pior ponto,

for inviável (



), volte para o passo 2 (ou modifique

para ser viável).

Calcule

(

). Se

é insatisfatório (



), volte para o passo 2.

Atualize o conjunto

substituindo o atual pior ponto pelo ponto tentativa: (





{

} / {

}). Encontre

na nova população

. Se

, então tome

como os novos

Se um critério de parada for satisfeito, pare; senão, volte para o passo 2.

As duas principais diferenças entre os CRSA disponíveis referem-se a (

) o modo

de geração do ponto tentativa (passo 2) e (

) o acesso opcional a uma fase de busca

local sempre que o melhor ponto for o mais recente na população (quando



passo 5). Deve-se notar que todas as versões assumem que



; como regra geral,

sugere-se tipicamente

= 10(



1).

5.3.2 – Algumas versões

O CRSA proposto por Price (1977) foi aparentemente o primeiro a apresentar o

formato acima. Ele não inclui uma fase local opcional no passo 5. O ponto tentativa no

passo 2 é gerado da seguinte maneira: escolha aleatoriamente



1 pontos distintos da

população atual

, ...,

(formando-se um

simplex



). Determine o centróide

dos

primeiros pontos

, ...,

. Em seguida, determine o ponto tentativa

como a

reflexão do ponto remanescente

em torno do centróide

= 2



(5.4)

Nesse esquema, a escolha do ponto

é puramente geométrica; ela não leva em

consideração o comportamento da função em torno de

. Apesar de sua objetividade,

boa varredura de busca e pequenos requisitos de memória, o CRSA de Price (1977)

exibe baixas taxas de convergência – o que motivou investigações posteriores por

melhorias.

Ali

et al.

(1997) compararam algumas versões de CRSA, numerando-as

cronologicamente como se segue: o CRS1, o algoritmo original de Price recém descrito

acima; o CRS2, também proposto por Price, fazendo-se um uso mais elaborado dos

simplexes

na geração do ponto tentativa; o CRS3, também devido a Price, no qual uma

fase local foi incluída; o CRS4, uma modificação do CRS2 incluindo-se buscas locais

aleatórias em torno do melhor ponto usando-se distribuição





; o CRS5, excluído das

comparações de Ali

et al

. (1997) porque apresenta uma busca local baseada em

gradientes e a ênfase fora dada a algoritmos que não usam gradientes; finalmente, o

CRS6, proposto por Ali

et al

. (1997), no qual as buscas locais do CRS4 (baseadas em

distribuição





) são mantidas e a fase global usa um esquema de interpolação

quadrática.

No esquema de interpolação quadrática do CRS6, três pontos distintos são

tomados da população atual

: o melhor ponto,

, e dois outros aleatórios,

Variando-se

de 1 a

, interpolações parabólicas são construídas usando-se as

correspondentes

-ésimas coordenadas desses três pontos,

. A coordenada

ponto tentativa é igualada ao ponto

extremo

(ponto crítico) da parábola. A idéia

heurística por trás desse esquema é considerar promissora qualquer região em torno do

melhor ponto e promover uma busca global por melhores pontos nessas regiões usando-

se uma interpolação de coordenadas sensata.

Um detalhe, porém, aparentemente não atraiu a atenção dos autores do CRS6.

No artigo de Ali

et al.

(1997), é dito que o componente

do ponto tentativa é o

mínimo

da parábola que interpola as

-ésimas coordenadas daqueles três pontos. Entretanto, essa

parábola não apresenta necessariamente um mínimo; ao contrário, o extremo pode ser

um máximo, Figura 23.

Parece que a expectativa de Ali

et al.

(1997) era de que apenas situações de

mínimo pudessem ocorrer. Daí, supondo-se que o comportamento da função objetivo

projetada no plano



não seja muito diferente de uma parábola, compreende-se o

porquê de se tomar o mínimo como uma boa escolha para a coordenada do ponto

tentativa. O que Ali

et al.

(1997) talvez não perceberam é que situações de máximo, em

geral, podem ser freqüentes nesse esquema (as probabilidades de ocorrência de

mínimos e máximos são normalmente da mesma ordem), contradizendo o argumento

dado ainda há pouco. Os máximos, entretanto, provavelmente tornam a busca mais

exploratória, carregando às vezes o ponto tentativa para longe do atual melhor ponto,

evitando-se assim convergências prematuras para mínimos locais. Portanto, o esquema

de interpolação usado no CRS6 possui basicamente dois aspectos: (

) quando o extremo

da parábola é um mínimo, a busca é mais local, mais intensificada, embora mesmo

quando isso ocorra, o ponto tentativa pode cair longe do atual melhor ponto (por

exemplo, se a interpolação for mal condicionada, o mínimo pode ficar muito afastado

dos três pontos usados na interpolação); (

) quando o extremo é um máximo, a busca é

mais global, mais diversificada. Esse aspecto não estava muito claro no artigo de Ali

al.

(1997) e acredita-se que essas duas características sejam ambas responsáveis pelo

bom desempenho do algoritmo.

Figura 23 – Interpolações quadráticas do CRS6/CRSI. Os extremos das parábolas

podem ser mínimos ou máximos.

Valor objetivo

Quadráticas

Projeção da função

Extremo - coordenada do ponto tentativa

As comparações de Ali

et al.

(1997) foram realizadas em 13 problemas-teste

retirados de referências clássicas sobre otimização global. O critério de parada adotado

foi







. Em todos os testes, os resultados mostraram um desempenho superior

do CRS6 sobre as outras versões de CRSA. Também concluiu-se que as interpolações

quadráticas da fase global foram mais propícias à melhoria de desempenho do que as

buscas locais baseadas em distribuição





(embora elas também tenham exercido um

efeito benéfico). Em outro artigo, Ali e Törn (2004) compararam os desempenhos de

alguns algoritmos populacionais para otimização global. Essas comparações incluíram

alguns algoritmos genéticos, de evolução diferencial e de busca aleatória controlada,

neste último caso, as versões CRS2, CRS4, CRS6 e CRSI, essa última sendo

equivalente ao CRS6 mas sem usar buscas locais. Seus resultados mostraram um bom

desempenho do CRS6 e CRSI em termos de chamadas da função objetivo e tempo de

CPU, embora o mesmo não possa ser dito em termos de confiabilidade em se encontrar

mínimos globais.

Manzanares Filho

et al.

(2005) propuseram um CRS6 modificado e o aplicaram

ao problema de projeto inverso de grades de aerofólios. Essa nova versão, denominada

CRS-VBR (

Controlled Random Search using Variability Based Reflections

– Reflexões

Baseadas em Variabilidade), faz um uso seletivo das interpolações quadráticas do

CRS6, leva em consideração a variabilidade da função em torno do melhor ponto e

elimina as buscas locais. Como no CRS6, três pontos distintos são tomados da

população atual

: o melhor ponto,

, e dois outros aleatórios,

. Seus respectivos

valores de função objetivo são representados por

e suas respectivas

-ésimas

coordenadas por

. Um valor médio,

, e uma medida de variabilidade local,



em torno do melhor ponto são calculados por:





(5.5)







(5.6)

onde

é o valor objetivo do atual pior ponto da população.

A geração do ponto tentativa proposta por Manzanares Filho

et al

. (2005) é

descrita como se segue:

Faça um teste de condicionamento da interpolação quadrática entre

Se ela for considerada mal condicionada, então calcule

aleatoriamente no

correspondente intervalo viável:



Senão, se a coordenada do melhor ponto situa-se entre as coordenadas dos

outros dois pontos – i.e., se (

–

)(

–

)



0 – então calcule

como o

mínimo da interpolação quadrática entre

, como no CRS6.

Senão, calcule um centróide ponderado cujas coordenadas

são dadas por:

)()(

3322

jljl

ffff

rffrff











(5.7)

e daí calcule a coordenada

do ponto tentativa usando-se uma reflexão

baseada em variabilidade de

em torno de

jjj

glp )1()2(









(5.8)

Nesse esquema, o cálculo aleatório de

quando a interpolação é considerada

mal condicionada incrementa a varredura da busca e ajuda a evitar uma contração

prematura da população em direção a um mínimo local. Até mesmo uma eventual falha

do processo numérico pode ser evitada com essa medida. No passo 2, satisfazer a

condição de aceitação para a coordenada do ponto tentativa do CRS6 realmente

assegura que essa coordenada corresponda ao mínimo da interpolação e que também

seja viável. O cálculo da coordenada

do centróide na Equação (5.7) coloca mais peso

na coordenada cujo valor objetivo é maior, de modo a enfatizar a variabilidade local em

torno do ponto mínimo.

O cálculo da coordenada do ponto tentativa na Equação (5.8) com uma medida

de variabilidade local (



) – Equações (5.5) e (5.6) – representa um procedimento

heurístico para automaticamente balancear as buscas global e local. Para baixa

variabilidade local (





0), a coordenada do centróide é refletida para longe da

coordenada do melhor ponto, pois o correspondente intervalo de interpolação não

aparenta ser promissor. Isso ajuda o algoritmo a escapar de pontos de mínimo local e

torna a busca “mais global” quando necessário. Para alta variabilidade local (





1) o

centróide é refletido próximo à coordenada do melhor ponto, pois melhores valores da

função objetivo agora são esperados em torno de

. Isso ajuda o algoritmo a acelerar a

busca por mínimos locais quando necessário.

O CRS-VBR mostrou-se mais eficiente do que as versões CRS6 ou CRSI no

problema de projeto inverso de grades de aerofólios (Manzanares Filho

et al

., 2005),

trazendo a expectativa de perfazer com êxito também o projeto otimizado de turbinas

hidráulicas axiais deste trabalho.

Um CRSI modificado também foi testado antes de ser aplicado ao problema da

turbina axial. Esta versão consiste no CRSI original de Ali

et al

. (1997) com duas

pequenas alterações: (

) um teste de mal condicionamento para as interpolações

quadráticas é feito

a priori

(no CRSI ou CRS6, o condicionamento das interpolações

não é checado); (

) quando a coordenada do extremo da parábola situa-se fora do

correspondente intervalo viável, a coordenada do ponto tentativa é igualada à fronteira

adjacente, uma forma de viabilização forçada no passo 3 do CRSA básico. Dessa forma,

a geração do ponto tentativa deste CRSI modificado é a seguinte:

Tome três pontos distintos da população atual

: o melhor ponto,

, e dois

outros aleatórios,

. Suas

-ésimas coordenadas são representadas

respectivamente por

, (

= 1 , ...,

) e os valores de função objetivo

correspondentes, por

Faça um teste de condicionamento da interpolação quadrática entre

Se for considerada mal condicionada, então calcule

aleatoriamente no

correspondente intervalo viável:

xpx 

Senão, calcule

como o

extremo

da parábola passando por

322332

)()()(

frlflrfrr

jjjjljj





(

= 1, ...,

) (5.9)

for inviável (



), porque algum



];[

, então:

(

) se



, faça

xp 

;

(

) senão, se



, faça

xp 

Observações



A interpolação pode ser considerada mal condicionada quando o módulo do

determinante associado ao termo quadrático da parábola,

322332

)()()( frlflrfrr

jjjjljj



, for inferior a uma pequena tolerância

. No programa em MatLab



, fez-se

eps

, onde

eps

é a precisão decimal

relativa do computador – distância entre 1,0 e o próximo maior número

decimal (em dupla precisão de 32 bits,

eps

= 2,2204







Como no CRS-VBR, quando a interpolação é considerada mal condicionada,

o cálculo aleatório de

incrementa a varredura da busca e ajuda a evitar

contrações prematuras da população. Até mesmo uma eventual falha do

processo numérico pode ser evitada com essa medida.



Fazer a coordenada

do ponto tentativa igual ao limite adjacente do

intervalo correspondente quando

é inviável (passo 4 da geração do ponto

tentativa) evita que se recalcule todas as outras coordenadas (como ocorre no

CRSI/CRS6 originais; passo 3 do CRSA básico). Além do que, isso ajuda a

levar a população para as fronteiras do “hiper cubo”

quando necessário.

Antes de se aplicar algum algoritmo de busca aleatória controlada ao problema

de projeto da turbina hidráulica axial, decidiu-se avaliar o desempenho de três versões

de CRSA em alguns problemas-teste da literatura sobre otimização global, Tabela 2

(Dixon e Szegö, 1978; Ali e Törn, 2004; Price, 1977). Os CRSA’s foram

implementados em linguagem de programação MatLab



, pois o programa de análise do

escoamento na turbina axial também foi escrito em MatLab



. As versões de CRSA

investigadas são o CRSI modificado (que acabou de se descrever), o CRS-VBR e uma

variante do CRS-VBR, denominada CRS-VBR-Explor, não descrita aqui mas incluída

nas comparações para ilustração (o CRS-VBR-Explor tenta melhorar o caráter global da

busca do CRS-VBR, usando-se mais aleatoriedade no começo do processo e mais

reflexões no fim). Cada caso foi rodado pelo menos trinta (30) vezes.

Na Tabela 2, listam-se os valores médios de chamadas da função objetivo, FE –

essa é uma medida de rapidez do algoritmo – e sucessos totais, TS, i.e., quando um

mínimo global é encontrado – essa é uma medida de confiabilidade ou robustez do

algoritmo. Para o CRSI, também são listados os valores médios das relações entre as

ocorrências de mínimos e máximos nas interpolações quadráticas (

min/max

). Note que a

idéia não é levantar estatísticas precisas, mas apenas evidenciar as características

principais desses algoritmos (Ali

et al

.(1997), por exemplo, consideraram suficiente

rodar apenas três (3) vezes cada um de seus problemas-teste para se obter um

comportamento médio confiável). Em todos os casos, usou-se

= 10(



1) e







como critério de parada.

Tabela 2 – Comparação de três versões de CRSA

Função Objetivo CRS-VBR CRS-VBR-Explor CRSI modificado

FE TS FE TS FE TS

min/max

Price [1] 300 10 % 500 40 % 430

70 %

Goldstein e Price (GP) [2] 220 90 % 220 90 % 160

95 %

Branin (BR) [2] 180 100 % 200 100 % 200 100 % 2

Shekel 5 (S5) [2] 400 40 % 450 30 % 600

45 %

Shekel 7 (S7) [2] 450 25 % 500 20 % 800

80 %

Shekel 10 (S10) [2] 400 10 % 470 10 % 900

90 %

Hartman 3 (H3) [2] 250 85 % 250 85 % 250

99 %

1.5

Hartman 6 (H6) [2] 550 50 % 700

80 %

900 70 % 1.3

Rastrigin (RG5) [3] Nenhum mínimo global foi encontrado

Rosenbrock (ER10) [3] Nenhum mínimo global foi encontrado

[1]: Price (1977)

[2]: Dixon e Szegö (1978)

[3]: Ali e Törn (2004)

O gráfico altamente “rugoso” da função bidimensional (Price), com 49 mínimos

locais de valor 1,0 e um único mínimo global, na origem, de valor 0,9, está

representado na Figura 24. Essa é uma das situações mais difíceis em otimização

global, pois o mínimo global é apenas ligeiramente diferente dos diversos e próximos

mínimos locais, sendo difícil de ser encontrado com alta taxa de sucesso. Nesse

problema, o melhor algoritmo foi o CRSI, com 70% de TS e 430 FE (o CRS1 original

de Price (1977), quando encontra o mínimo global, o faz em média com 700 FE). Um

algoritmo de busca local, como o SQP do MatLab





fmincon

fatalmente convergiria

para algum dos muitos mínimos locais, a menos que o ponto de partida

já estivesse

suficientemente próximo do mínimo global (origem).

Deve-se notar que excetuando-se (Price), (RG5) e (ER10), os outros problemas

são considerados como fáceis por Ali e Törn (2004). No entanto, mesmo para alguns

desses problemas fáceis, os algoritmos testados encontraram dificuldades em encontrar

um mínimo global (S5, S7, S10, H6). Os problemas moderadamente difíceis (RG5) e

(ER10) simplesmente não foram resolvidos por nenhum dos três CRSA’s.

A Tabela 2,

porém, aponta o CRSI como o mais exploratório e eficiente entre eles, apresentando as

maiores taxas de sucesso e um número razoável de chamadas de função objetivo

. Outra

observação é que os extremos das interpolações quadráticas do CRSI freqüentemente

são máximos em vez de apenas mínimos, como se sugeriu em Ali

et al.

(1997); em

geral, seus números de ocorrência são da mesma ordem. Note também que o valor

min/max

diminui com o grau de dificuldade dos problemas, medido pelos sucessos

totais TS.

Figura 24 – Função bidimensional de Price (1977). O mínimo global, levemente

diferente dos diversos mínimos locais, está na origem. Problemas

multimodais estão entre os mais difíceis em otimização global.

Dessa comparações, ponderou-se o CRSI modificado como o mais promissor

entre as versões testadas para executar a otimização do projeto de turbinas hidráulicas

axiais, Figura 25. De fato, alguns testes preliminares foram realizados com esses

algoritmos no problema da turbina axial, e o CRSI modificado mostrou-se realmente o

mais eficiente entre eles, convergindo com maior freqüência para a solução global (que

ainda será apresentada) e o fazendo com menos chamadas de função objetivo



programa de análise do escoamento, desenvolvido no capítulo 4. Esse último aspecto é

especialmente relevante na otimização de turbomáquinas, onde são comuns

solvers

alto custo computacional.

Figura 25 – Fluxograma do projeto otimizado usando-se CRSA.

Restrições laterais, x e x (região viável S),

população inicial

P de N pontos aleatórios em S

Otimizador (CRSA)

Ciclo interno heurístico (CRSI)

em direção a um ótimo da função objetivo

Função objetivo

solver

Critério de parada

satisfeito?

Projeto otimizado

Variáveis de projeto,

Rendimento, 

Altura, H

NÃO

SIM

Capítulo 6

RESULTADOS E DISCUSSÃO

6.1 ANÁLISE DO PROJETO INICIAL DE SOUZA (1989)

Em 1989, o Prof. Zulcy de Souza e sua equipe do LHPCH – UNIFEI projetaram,

construíram e ensaiaram um modelo de turbina hélice do tipo tubular, semelhante às

Figuras 2 e 3, para a antiga Mecânica Pesada de Taubaté – SP. Como estão disponíveis

as planilhas de projeto e os resultados de ensaio dessa turbina (Souza, 1989), ela foi

escolhida para servir de comparação com os resultados de otimização, sendo referida

aqui como o “projeto inicial”.

A Tabela 3 apresenta algumas características básicas desta turbina, cuja seção

meridional é esquematicamente a da Figura 3, repetida a seguir como Figura 26 para

maior comodidade. Seu projeto foi feito segundo a hipótese de “vórtice potencial”

(discutida no item 4.1) e usando-se a Teoria da Asa de Sustentação (ou Teoria do

Elemento de Pá) para se determinar os perfis das pás e seus ângulos de montagem –

uma metodologia de projeto clássica (Bran e Souza, 1979). Foram escolhidos perfis

Göttingen (tradicionais em turbinas hidráulicas) em 5 estações radiais equiespaçadas,

sem se fazer uso de qualquer tipo de parametrização. No entanto, as parametrizações

parabólicas adotadas neste trabalho reproduzem satisfatoriamente a geometria do

projeto inicial, conforme visto no item 3.2, Figura 9.

Antes, porém, de executar as rotinas de otimização, o projeto inicial de Souza

(1989) foi analisado usando-se a metodologia de cálculo do escoamento desenvolvida

no capítulo 4. Essa metodologia foi implementada num programa em MatLab

acordo com o fluxograma geral da Figura 21 e demais explicações dadas no capítulo 4.

Tabela 3 – Características básicas do projeto inicial (Souza, 1989)

Vazão, Q

0,286 m³/s

Rotação, n

1145 rpm

Altura, H

4,0 m

Rotação específica,

4/3

2/1

)(gH

652

Rendimento esperado,

85 %

Potência de eixo esperada, P

9,5 kW

Diâmetro externo, D

280 mm

Diâmetro interno, D

112 mm

Número de pás, N

pá

Ângulo de saída do distribuidor,

60°

Figura 26 – Esquema da seção meridional da turbina hélice projetada por Souza (1989).

Distribuidor Rotor

Tubo de sucção

Inicialmente, os fatores empíricos presentes na formulação de perdas (Tabela 1)

foram ajustados de modo a se reproduzir razoavelmente as características de

funcionamento da turbina ensaiada. Não se trata de uma validação no sentido exato do

termo, mas apenas de um ajuste particular para a turbina do projeto inicial. A Tabela 4

apresenta então os valores adotados para esses fatores empíricos. Na Tabela 5,

comparam-se os dados de projeto, os resultados ótimos do ensaio (Souza, 1989) e

aqueles fornecidos pelo programa de análise (com os valores da Tabela 4 para os

fatores empíricos). Deve-se ter em mente que o objetivo aqui não é reproduzir

exatamente as características de funcionamento medidas, mas apenas obter um simples

ajuste que indique as corretas

tendências

de perdas, pois esse é um aspecto fundamental

para o êxito do processo de otimização.

Tabela 4 – Valores adotados para os fatores empíricos de perdas

Símbolo Descrição

Valor Escolhido

[

Faixa

]

Coeficiente de perda por choque na entrada do rotor

0,4 [0,5 - 0,7]

Coeficiente de perda de perfil em grades motoras para

deflexão nula, na correlação de Soderberg

0,05 [0,04 - 0,06]

Fator de perda por difusão no tubo do sucção 0,09 [0,09 ou 0,12]

Fator de perda por giro no tubo de sucção 0,4 [0,2 - 0,4]

mec

Rendimento mecânico ou externo 0,96 [0,96 - 0,99]

Tabela 5 – Dados de projeto, resultados ótimos do ensaio e resultados do programa

para o projeto inicial de Souza (1989)

Par vazão – rotação

Grandezas Resultantes

Condições

Projeto

Condições Ótimas

Ensaio

Resultados

Programa

Vazão,

/s) 0,286 0,267 0,288

Rotação,

(rpm) 1145 1145 1145

Altura,

(m) 4,0 3,9 3,7

Rendimento,

(%)

85,0 82,0 80,7

Potência de eixo,

(kW) 9,5 8,3 8,4

O par vazão – rotação escolhido no programa,

= 0,288 m³/s e

= 1145 rpm,

não é o de projeto nem o ótimo de ensaio (Tabela 5). Porém, com a geometria do

projeto inicial de Souza (1989), produz-se

= 3,7 m,

= 80,7 % e

= 8,4 kW,

valores esses considerados moderadamente próximos das medições. De fato, seria mais

correto proceder a um ajuste mais completo dos resultados do programa com as

medições, cobrindo uma certa faixa de operação da turbina. Para os propósitos de

agora, no entanto, isso não chega a ser indispensável e acredita-se que o ponto

–

escolhido, juntamente com os valores da Tabela 4 para os fatores empíricos, sejam

razoáveis para que as comparações entre a geometria otimizada e o projeto inicial

tenham sentido.

Com essas considerações em mente, apresenta-se a Tabela 6, com resultados

mais completos da análise do projeto inicial. Em seguida, analisam-se as Figuras 27 a

30, as quais trazem as distribuições radiais de várias grandezas importantes do

escoamento associadas ao campo de velocidade, determinadas segundo os

desenvolvimentos do capítulo 4.

Tabela 6 – Análise do projeto inicial

Variáveis de Projeto

Grandezas Resultantes

Projeto inicial (Souza, 1989)

com a geometria parametrizada

cubo meio ponta

(°)

49,3 26,2 17,4

(

−

)

1,61 1,08 0,889

(%)

considerou-se pás retas

(°)

60,0

pá

(W) 8761

(W) 215

(

r + ts

)

(W) 1443

(%)

84,09

(%)

80,72

(m) 3,70

A Figura 27 apresenta os perfis de velocidade nas seções de saída do distribuidor

e do rotor, calculadas de acordo com a condição de equilíbrio radial (itens 4.2 e 4.3).

Lembre-se que essa cinemática é a base de todos os cálculos subseqüentes para a

avaliação de desempenho da turbina. Representa-se também, em linhas tracejadas, a

hipótese de vórtice potencial usada por Souza (1989), segundo a qual as velocidades

meridionais (

) seriam uniformes, o giro na saída do distribuidor (

) teria

variação hiperbólica (

= const.

) e o giro na saída do rotor (

) seria nulo. Como se

percebe, essas hipóteses de projeto não se verificaram na prática, principalmente porque

o distribuidor dessa turbina, tendo um ângulo de saída fixo ao longo do raio, não

propicia o vórtice potencial no recinto entre distribuidor e rotor. Isso reforça a

importância de se considerar corretamente o equilíbrio radial para uma predição

fisicamente mais consistente dos perfis de velocidade.

O rotor do projeto inicial deflete o escoamento de tal modo que a distribuição de

fica mais variável que a de

) e o giro na saída (

) assume valores

ligeiramente positivos ao longo de todo o raio (

90°), mas não a ponto de causar

perdas excessivas (por giro) no tubo de sucção.

Figura 27 – Perfis de velocidade no Figura 28 – Distribuições de

projeto inicial.

Na Figura 28, representam-se as distribuições de

, a quantidade de movimento

angular (axial) por unidade de massa do fluido. Lembrar que a diferença entre

) e

relaciona-se diretamente com o trabalho específico das pás (Equação de

Euler das turbomáquinas, (4.43)) e, conseqüentemente, tem papel fundamental no

desempenho da turbina. Novamente, fica claro que a hipótese de vórtice potencial (

0.02.04.06.08.0

Velocidade (m/s)

0.4

0.6

0.8

1.0

r/rt

0.02.04.06.08.0

Vórtice

potencial

-0.10 0.00 0.10 0.20 0.30 0.40

Momento cinemático (m²/s)

0.4

0.6

0.8

1.0

r/rt

-0.10 0.00 0.10 0.20 0.30 0.40

Vórtice

potencial

= const

) não foi verificada, como já se explicou (distribuidor sem torção radial).

Mencione-se também que a Figura 28 mostra os valores calculados nas

(= 7) grades

escolhidas do cubo à ponta das pás e que para representar a variação de

(necessária

à Equação (4.37)) escolheu-se uma regressão parabólica (item 4.2) e para

necessária na Equação (4.53), uma cúbica (item 4.3). Note mais uma vez que o giro na

saída do rotor é positivo ao longo de toda a extensão radial.

A Figura 29 traz as distribuições das perdas específicas que ocorrem na turbina

original, segundo a modelagem adotada no item 4.5. Como é de se esperar numa turbina

hidráulica, as perdas no distribuidor (onde o escoamento é bastante acelerado),

, são

as menos acentuadas. Além disso, por serem pouco variáveis, justifica-se sua

desconsideração na equação de equilíbrio radial após o estator (Equação (4.34)). As

perdas de perfil no rotor,

, crescem do cubo para a ponta essencialmente porque, em

turbinas hidráulicas axiais, a corda

dos perfis também aumenta bastante nesse sentido

(veja a Figura 4), ou seja, há de fato mais atrito na região da ponta do que na região da

raiz das pás e esse efeito é considerado na correlação de Soderberg.

Devido à forma simplificada como se tratou a questão do choque na entrada do

rotor, fazendo-se o componente de choque igual ao componente de incidência (Equação

(4.72)), as perdas por choque,

Pch

, ficaram demasiado elevadas no projeto inicial,

especialmente na região da ponta, motivo pelo qual adotou-se para o coeficiente de

choque

um valor pouco abaixo da faixa recomendada (Tabela 4). De fato, a turbina

original de Souza (1989) não deve operar em condições de choque tão exageradas,

como sugere a Figura 29. Porém, como se frisou no item 4.5, o importante é que,

mesmo sem uma predição quantitativamente acurada, esse modelo indique para o

otimizador que, com níveis menores de incidência, é possível reduzir-se as perdas por

choque, e essa é uma

tendência

correta.

As perdas no tubo de sucção,

Pts

, aparentam um comportamento razoável,

respondendo por cerca de metade das perdas totais, como realmente se espera em

turbinas hidráulicas de baixa queda (Drtina e Sallaberger, 1999; Moura e Brasil Júnior,

2005). Isso mostra que o modelo de Raabe (1985) é satisfatório para nossos propósitos

de otimização e que o ajuste para os fatores

(Tabela 4) está apropriado.

Observar que o nível superior de perdas na região do cubo deve-se às grandes

velocidades meridionais na saída do rotor (

) nessa região (Figura 27 e Equação

(4.74)).

A soma das perdas no rotor e no tubo de sucção,

(

)

também está

representada na Figura 29, pois essas perdas são integradas no mesmo fluxo de massa

(saída do rotor). Mencione-se que a Figura 29, como a Figura 28, também apresenta os

valores calculados em

(= 7) estações radiais. No item 4.5.6, escolheram-se

polinômios do terceiro grau para representar as variações de

(

)

, pois conferem

um bom ajuste a esses perfis.

Figura 29 – Perdas no projeto inicial. Figura 30 – Trabalho específico no

projeto inicial.

Por fim, a Figura 30 traz a distribuição de trabalho específico nas pás do rotor,

pá

, onde se percebe mais uma vez que a condição de vórtice potencial (

pá

= const.)

não foi cumprida. Percebe-se também um carregamento bem maior na ponta do que na

raiz das pás (

pát

páh

≅

1,8 ), o que pode inclusive ser perigoso para a cavitação. Esse

comportamento, aliás, depreende-se diretamente da Figura 28, onde as curvas de

afastam-se à medida que se caminha do cubo para a ponta.

Os resultados globais finais (integrados) de perdas, potências e rendimentos para

o projeto original já foram antecipados na Tabela 6.

De toda essa análise, pode-se concluir que, embora as hipóteses de projeto

(vórtice potencial) não tenham sido exatamente satisfeitas pelo campo de escoamento

resultante da geometria original de Souza (1989), esse projeto inicial já é bem razoável

e não deve estar muito afastado de uma configuração ótima.

0.0 2.0 4.0 6.0

Perda específica (J/kg)

0.4

0.6

0.8

1.0

r/rt

0.0 2.0 4.0 6.0

Pch

Pts

P(r + ts)

0 10203040

Trabalho específico (J/kg)

0.4

0.6

0.8

1.0

r/rt

0 10203040

Vórtice Potencial

6.2 RESULTADOS DO SQP-FMINCON

Ao menos a princípio, poder-se-ia executar uma rotina de otimização definindo-

se somente o par vazão – rotação e as restrições de altura, e deixando-se a geometria da

turbina inteira para ser tratada como variáveis de projeto. Além de existir, muito

provavelmente, uma grande diversidade de mínimos locais para esse problema, analisar

uma solução desse tipo seria bastante difícil, pois estariam sendo tratadas

simultaneamente dimensões básicas da turbina, como diâmetros do rotor e número de

pás, e dimensões já mais detalhadas, como ângulos de montagem e arqueamento dos

perfis.

De fato, para o nível de sofisticação proposto neste trabalho, otimizar toda a

geometria da máquina é mesmo impraticável, pois a modelagem desenvolvida não é

capaz de considerar corretamente todos os efeitos que estariam envolvidos.

Muito mais construtivo aqui é preestabelecer as dimensões “macroscópicas”

principais da turbina e deixar como variáveis de projeto apenas aquelas geometrias

“microscópicas” tratadas no capítulo 3, até porque, para as dimensões mais básicas, já

existem diretrizes de projeto realmente otimizadas e que podem ser usadas no início do

dimensionamento (Schweiger e Gregori, 1988).

Resolveu-se, então, usar aquela mesma turbina projetada por Souza (1989),

mantendo-se toda a sua geometria, à exceção das dimensões tratadas como variáveis de

projeto. Dessa forma, pode-se comparar diretamente as soluções fornecidas pelo

programa de otimização com o projeto inicial. Isso vai possibilitar, inclusive, uma

melhor crítica para a modelagem e o procedimento de cálculo adotados.

Sendo assim, as otimizações usando-se o SQP-

fmincon

foram realizadas de

acordo com a Tabela 7, para o ponto de projeto e restrições de altura, e Tabela 8, para

as restrições de faixa das variáveis de projeto.

Tabela 7 – Ponto de projeto e restrições de altura para a otimização

Vazão volumétrica, Q

0,288 m³/s

Rotação, n

1145 rpm

Altura mínima, H

3,5 m

Altura máxima, H

3,7 m

Tabela 8 – Restrições de faixa para as variáveis de projeto

Variável de

projeto

(°)

cubo meio ponta

(

−

)

cubo meio ponta

(%)

cubo meio ponta

(°)

40 25 15 1,61 1,08 0,889 0,8 0,5 0,1 50

55 35 25 1,70 1,20 1,00 6,0 4,0 2,0 70

A Tabela 9 apresenta três soluções encontradas pelo

fmincon

, partindo-se de

diferentes estimativas iniciais

, juntamente com seus respectivos parâmetros de

desempenho. Por se tratar de um método baseado em gradiente, o critério de parada

nesses casos é uma tolerância para a derivada direcional: aceita-se como mínimo local

um ponto onde a direção de busca – derivada direcional – tenha norma inferior a

−

. Esse valor depende da função objetivo em questão e deve ser ajustado de modo

que não haja uma interrupção prematura do processo de convergência e nem cálculos

em demasia da função objetivo. As soluções apresentadas exigiram de 133 a 189

chamadas do

solver

(função objetivo), valores esses considerados razoáveis e

praticáveis de se perfazer em um único PC, dentro do espírito de metodologias de baixo

custo computacional explicadas nos capítulos 1 e 3. Por exemplo, rodando em um

Pentium4

™

3,0 GH

, essas otimizações levaram cerca de apenas dois (2) minutos.

Tabela 9 – Três soluções encontradas usando-se o SQP-

fmincon

Variáveis de Projeto

Grandezas Resultantes

fmincon

= (49,3 26,2 17,4

1,61 1,08 0,889

0 0 0 % 60)

cubo meio ponta

fmincon

= (55,0 32,0 25,0

1,61 1,08 0,889

0,8 0,5 0,2 % 51)

cubo meio ponta

fmincon

= (40,0 25,0 15,0

1,61 1,08 0,889

0,8 0,5 0,2 % 68)

cubo meio ponta

(°)

49,4 28,7 18,7 55,0 30,1 23,3 40,5 27,4 17,5

(

−

)

1,61* 1,08* 0,889* 1,614 1,08* 0,889* 1,61* 1,08* 0,889*

(%)

3,64 1,67 0,93 5,87 2,17 0,1* 3,34 1,22 1,19

(°)

60,7 50,5 68,9

pá

(W) 9218 9046 9208

(W) 211 278 178

(

)

(W) 1011 1040 1104

(%)

84,77 83,79 84,27

(m) 3,70 * 3,70 * 3,70 *

(*) Observação: os valores marcados com um asterisco na Tabela 9 correspondem a uma

restrição que foi ativada na solução final. Veja as Tabelas 7 e 8.

A primeira solução – ‘

fmincon

1’ – foi gerada partindo-se do projeto original

para as variáveis de projeto (Tabela 6), projeto inicial esse já bastante razoável,

conforme se analisou no item 6.1. Para as demais soluções, partiu-se de pontos

diferentes, mas que também são minimamente favoráveis ao equilíbrio radial. De fato,

geometrias desfavoráveis ao equilíbrio radial não podem ser analisadas usando-se o

programa desenvolvido, pois nesses casos não se consegue a convergência dos perfis de

velocidade na saída do rotor. Isso se deve à forma como foi imposta a condição de

equilíbrio radial na saída das grades (itens 4.2 e 4.3). Tal fato, todavia, não deve ser

encarado como uma limitação séria da metodologia, afinal o objetivo é projetar turbinas

otimizadas, nas quais, normalmente, o estabelecimento do equilíbrio radial é favorecido

pela geometria sem maiores dificuldades.

Deve-se notar os distintos valores para o ângulo de saída ótimo do distribuidor

(

) entre as três soluções. Essa é uma questão central porque o escoamento na saída do

distribuidor afeta fortemente a geometria restante das pás do rotor.

De fato, a solução ‘

fmincon

1’ é o melhor projeto que foi encontrado partindo-

se o otimizador SQP-

fmincon

de vários pontos iniciais

diferentes. Isso sugere a

possibilidade de ‘

fmincon

1’ ser o ótimo global de nosso problema.

A segunda solução

−

‘

fmincon

2’

−

possui um ângulo de saída do distribuidor

muito pequeno (

= 50,5°). Assim, esse distribuidor confere muita quantidade de

movimento angular para o escoamento na entrada do rotor. Devido ao limite superior

para a altura de energia (

= 3,7 m), o rotor não pode extrair toda essa quantidade de

movimento angular e o fluxo na saída do rotor dá-se com um giro positivo significativo

(

90°; i.e., giro no sentido de rotação das pás), causando menos absorção de

potência nas pás e mais perdas no rotor e no tubo de sucção em comparação às da

solução ‘

fmincon

1’ (Tabela 9). Além disso, a deflexão excessiva nas palhetas

diretrizes aumenta consideravelmente as perdas no distribuidor (32 %). Como resultado

final, o rendimento do projeto ‘

fmincon

2’ é 1 % inferior que o de ‘

fmincon

1’.

A solução ‘

fmincon

3’, por outro lado, possui um ângulo de saída do

distribuidor muito grande (

= 68,9°). Nesse caso, o escoamento na entrada do rotor

apresenta pouca quantidade de movimento angular. Daí, para que se produza uma

potência de eixo razoável, as pás do rotor defletem o escoamento de modo que o giro na

saída torna-se bastante negativo (

90°), aumentando-se em 10 % as perdas no rotor

e no tubo de sucção quando comparadas ao projeto ‘

fmincon

1’. Embora as perdas no

distribuidor agora sejam menores, o rendimento é 0,5 % inferior ao de ‘

fmincon

1’.

A solução ‘

fmincon

1’ para o distribuidor (

= 60,7°) provavelmente leva à

quantidade de movimento angular ótima para a entrada do rotor. Isso torna possível a

absorção da energia do fluido pelas pás do rotor da maneira mais eficiente. De fato, os

ângulos de montagem e os arqueamentos dos perfis devem ser tais que se tenha choque

mínimo na entrada do rotor, assim como deflexões razoáveis juntamente com giro

mínimo na saída (

≈

90°). A modelagem adotada para as perdas e os desvios, com a

parametrização da geometria e a imposição da condição de equilíbrio radial, realmente

permitem ao otimizador seguir essas tendências corretas de projeto de turbinas

hidráulicas. O projeto inicial de Souza (1989), ponto

de ‘

fmincon

1’, constitui uma

alternativa razoável principalmente devido à escolha muito boa para

(= 60°). Essa

turbina, porém, ainda é susceptível de melhorias na geometria do rotor.

(a) (b)

Figura 31 – Ângulo de montagem (a) e arqueamento (b) nos projetos inicial e ótimo.

Como será mostrado no próximo item, ‘

fmincon

1’ é muito provavelmente o

ótimo global desse problema e por isso será referido daqui por diante como o “projeto

otimizado”. Na Figura 31, comparam-se os rotores original e otimizado. A principal

diferença está no ângulo de montagem da corda,

, que foi levemente aumentado ao

0.0 0.2 0.4 0.6 0.8 1.0

Raio adimensional (r/ rt)

Ângulo de montagem da corda, graus

Projeto otimizado

Projeto inicial

0.0 0.2 0.4 0.6 0.8 1.0

Raio adimensional (r/ rt)

-0.01

0.00

0.01

0.02

0.03

0.04

0.05

Razão flecha / corda (ARC)

Projeto otimizado

Pá reta

-0.01

0.00

0.01

0.02

0.03

0.04

0.05

longo de toda a extensão radial

−

as pás assim ficando mais abertas. O arqueamento

também é ajustado pelo otimizador de modo a se produzir pouca incidência (pouco

“choque”) na entrada do rotor (Figura 35), pouco giro na saída (Figuras 33 e 34) e bons

valores de trabalho específico (Figura 36). A pá do rotor otimizado, com seus ângulos

de montagem e arqueamentos, está representada na Figura 32.

Figura 32 – Perfis otimizados para as pás do rotor; as espessuras são apenas ilustrativas.

Figura 33 – Perfis de velocidade. Figura 34 – Perfis de momento angular por

unidade de massa do fluido.

cubo

meio

ponta

-0.10 0.00 0.10 0.20 0.30 0.40

Momento cinemático (m²/s)

0.4

0.6

0.8

1.0

r/rt

-0.10 0.00 0.10 0.20 0.30 0.40

Projeto inicial

Projeto otimizado

0.0 2.0 4.0 6.0 8.0

Velocidade (m/s)

0.4

0.6

0.8

1.0

r/rt

0.0 2.0 4.0 6.0 8.0

Projeto

inicial

otimizado

Nas Figuras 33 e 34 comparam-se, respectivamente, as distribuições de

velocidade e de quantidade de movimento angular específica (

) de ambos os projetos

(inicial e ótimo). Devido à diferença muito pequena entre os ângulos de saída do

distribuidor, os perfis de

são muito próximos entre os projetos inicial e

ótimo. Na saída do rotor, porém, os componentes de giro e meridional,

, são

diferentes. Enquanto a geometria original leva a uma grande variação radial de

, a

geometria otimizada produz um fluxo de saída mais uniforme, reduzindo-se a parcela

de difusão da perda no tubo de sucção. Outra observação importante é que, com a

geometria otimizada, o giro na saída do rotor (

) é negativo na região do cubo

e positivo na região da ponta, em vez de ser apenas positivo como ocorre no projeto

inicial. Essa tendência está em boa concordância com medidas de velocidade em

turbinas hidráulicas axiais bem projetadas (Osterwalder, 1960; Vivier, 1966; Lipej,

2004).

Figura 35 – Variação radial das perdas hidráulicas.

A Figura 35 apresenta a variação das perdas hidráulicas. As principais melhorias

estão nas perdas por choque (

Pch

) e no tubo de sucção (

Pts

). Como se percebe, o

otimizador buscou uma condição de incidência essencialmente nula na entrada do rotor.

Ainda que essa não seja exatamente a condição sem choque, fica claro que o conjunto

de variáveis de projeto escolhido dá liberdade suficiente para o otimizador alterar a

geometria buscando uma condição de escoamento favorável tanto na entrada como na

0.02.04.06.0

Perda específica (J/kg)

0.4

0.6

0.8

1.0

r/rt

0.02.04.06.0

Projeto

inicial

otimizado

P(r + ts)

Pch

Pts

saída do rotor. Embora a perda de perfil no rotor seja até maior no projeto ótimo, a

redução na perda por choque e no tubo de sucção é tal que

(

)

deve ser realmente

inferior no projeto otimizado, mesmo que

(

)

não alcance os 30 % de redução em

relação ao projeto inicial sugeridos pelas Tabelas 6 e 9 ou Tabela 11. De qualquer

maneira, essa análise sugere que, com modelos um pouco mais precisos, o processo de

otimização continuará seguindo as tendências de projeto corretas.

Na Figura 36, vê-se que o trabalho específico

pá

na solução ‘

fmincon

1’ é

superior e mais uniforme que o do projeto inicial. A maior uniformidade significa que o

carregamento hidrodinâmico é melhor distribuído nas pás otimizadas, algo importante

também para a questão da cavitação, não tratada aqui. Os maiores valores de

pá

foram

possíveis ao longo de todo o

span

porque as perdas são inferiores no projeto ótimo. Daí,

para um mesmo limite de altura de energia

disponível à turbina, o rotor pode extrair

mais potência do fluido (compare

pá

nas Tabelas 6 e 9; aumento de 5 %).

Figura 36 – Distribuição radial de trabalho específico.

Ainda que o ganho final de rendimento muito provavelmente não seja da ordem

de 4 % (Tabelas 6 e 9), o que para turbinas hidráulicas seria um valor excelente, é de se

esperar que o projeto otimizado apresente de fato um bom desempenho na prática,

mesmo que não possamos precisar qual o ganho de eficiência em relação ao projeto

inicial. O aspecto mais importante aqui é que a modelagem adotada e o procedimento

0 10203040

Trabalho específico (J/kg)

0.4

0.6

0.8

1.0

r/rt

0 10203040

Projeto inicial

Projeto otimizado

de cálculo desenvolvido, com a imposição do equilíbrio radial correto, permitam uma

boa

classificação

entre as soluções, pondo-se em evidência as tendências corretas para a

geometria e o campo de escoamento.

Uma dúvida natural que resta neste item é se há melhores projetos que o

‘

fmincon

1’, até aqui, a melhor solução apresentada (Tabela 9). De fato, partindo-se o

SQP-

fmincon

de alguma outra estimativa inicial

talvez seja possível alcançar

rendimentos ainda melhores. No próximo item, será apresentada a solução encontrada

usando-se o algoritmo de otimização global CRSI, indicando-se que ‘

fmincon

1’ seja

muito provavelmente o ótimo global desse problema.

6.3 RESULTADOS DO CRSA-CRSI

Aplicou-se o otimizador CRSI, descrito no capítulo 5, usando-se os mesmos

valores das Tabelas 7 e 8 para o ponto

–

de projeto, restrições de altura e região

viável das variáveis de projeto, afinal a idéia aqui é promover uma busca global no

mesmo espaço solução usado pelo SQP-

fmincon

As restrições para a altura de queda, como explicado no capítulo 2, são agora

impostas por meio de um esquema de penalização na função objetivo (Equação (2.4)).

Os estudos preliminares indicaram o valor

= 50 como adequado para o fator de

penalização. Essa escolha, como já foi explicado, não deve conduzir o processo de

otimização em direção a uma minimização apenas da penalidade, perdendo-se a

informação importante da função objetivo, i.e., o rendimento

. Por outro lado, as

restrições também não devem ser violadas ao fim do processo de convergência.

= 50

atende bem a essas exigências.

Outros dois parâmetros importantes do CRSA são o tamanho

da população e o

critério de parada. Para este último, usou-se a tolerância

−

, considerada

satisfatória para se evitar interrupções prematuras do processo de convergência assim

como um número exagerado de chamadas da função objetivo. Para o tamanho da

população, o valor genericamente sugerido

= 10(

1) = 10(10

1) = 110 mostrou-

se desnecessariamente grande, não trazendo benefícios adicionais em termos de

rendimento, mas apenas retardando a contração da população. Com apenas metade

desse valor (

= 55), obtém-se as mesmas soluções finais com um número até três (3)

vezes menor de chamadas do

solver

(função objetivo).

Devido ao caráter heurístico dos CRSA’s, qualquer solução deve ser tomada

segundo um senso estatístico: é preciso executar várias vezes a rotina de otimização

antes de se aceitar uma solução. Se o programa de análise do escoamento exigisse

grande esforço computacional e o método de otimização global chamasse muitas vezes

a função objetivo, seria proibitivo realizar um estudo estatístico ou mesmo executar

uma única rodada de otimização até que se atingisse o critério de parada. Rodando em

um Pentium4

™

3,0 GH

, essas otimizações gastam cerca de quinze (15) minutos, um

tempo relativamente pequeno e ainda dentro do espírito de metodologias de baixo custo

computacional.

Antes de apresentar os resultados de otimização global, deve-se fazer uma

consideração. Como se comentou no item 6.2, geometrias desfavoráveis ao equilíbrio

radial não podem ser analisadas usando-se o programa desenvolvido, já que nesses

casos não se consegue a convergência dos perfis de velocidade na saída do rotor. Por

causa disso, inclusive, para se inicializar o

fmincon

, a estimativa inicial

tinha de ser

minimamente favorável ao estabelecimento do equilíbrio radial na saída das grades. No

CRSA, porém, a população inicial de turbinas é tomada

aleatoriamente

na região viável

do problema. Daí fatalmente muitas dessas geometrias são bastante desfavoráveis ao

equilíbrio radial, algumas das quais produzindo refluxo e outras bombeando em vez de

“turbinar”, nenhuma delas podendo ser analisadas pelo programa. A maneira adotada

para contornar essa dificuldade foi a seguinte: se o número de

loops

do ciclo iterativo

para o perfil de velocidade na saída do rotor (Figura 18) chegar a 60, interrompe-se os

cálculos e atribui-se valor zero para o rendimento (

) e para a altura (

). Com esse

artifício, condenam-se as geometrias desfavoráveis ao equilíbrio radial – que também

não teriam bom desempenho –, pois o valor objetivo dessas geometrias torna-se grande

(devido à penalização) e elas acabam sendo inevitavelmente substituídas por qualquer

configuração mais favorável ao equilíbrio radial

−

com menos de 60

loops

no cálculo da

velocidade na saída do rotor. Mais uma vez, tal fato não deve ser encarado como uma

limitação séria da metodologia, afinal o objetivo é projetar turbinas otimizadas, nas

quais, normalmente, o estabelecimento do equilíbrio radial é favorecido pela geometria

sem maiores dificuldades.

A Tabela 10 apresenta três soluções obtidas usando-se o CRSA-CRSI. Observe

como elas são muito próximas entre si. De fato, outras soluções bastante semelhantes a

essas foram encontradas. A Tabela 10 também lista: o número de chamadas da função

objetivo (

), o número de iterações (i.e., número de vezes em que a população foi

atualizada,

iter

), o número de pontos tentativa gerados (

trial-p

), o número de vezes em

que as interpolações quadráticas foram mal condicionadas (

ill-cond

), o número de vezes

em que os extremos das parábolas foram mínimos (

min

) ou máximos (

max

), e o número

de vezes em que viabilizações forçadas foram realizadas (

viab

Tabela 10 – Soluções usando-se o CRSI modificado

Variável de Projeto

CRSI - 1 CRSI - 2 CRSI - 3

cubo

(°)

meio

ponta

49,0

29,7

18,7

48,7

28,7

19,9

49,5

29,0

18,6

cubo

(

−

)

meio

ponta

1,64

1,08 *

0,898

1,64

1,08 *

0,919

1,65

1,09

0,890

cubo

(%)

meio

ponta

3,22

1,94

0,85

3,86

1,66

1,20

3,69

1,77

1,02

(°)

58,7 59,6 60,5

(%)

84,66 84,61 84,71

(m) 3,70 * 3,70 * 3,70 *

921 1149 868

iter

553 699 493

trial-p

866 1094 813

ill-cond

4 9 0

min

4384 5154 3794

max

4204 5702 4292

viab

148 162 225

(*) restrição ativada; veja as Tabelas 7 e 8.

Deve-se notar os níveis de rendimento extremamente próximos alcançados pelas

três soluções. O limite superior para a altura de queda (

= 3,7 m) também sempre foi

ativado, tal qual ocorreu nos resultados do SQP-

fmincon

– para alturas maiores,

melhores rendimentos ainda podem ser obtidos. Um ângulo de saída do distribuidor

(

) de cerca de 60° parece ser mesmo o ótimo. Lembre-se que na melhor solução do

SQP-

fmincon

(

‘

fmincon

1’, Tabela 9), obteve-se

= 60,7°. A geometria do rotor

também está nitidamente determinada na Tabela 10; há apenas pequenas diferenças

entre as soluções CRSI para

, permitidas pelo

solver

do escoamento para o

mesmo nível de rendimento.

Alguns aspectos da otimização global merecem ser discutidos. Primeiro, o

número de chamadas da função objetivo, em torno de 1000, é razoável, sendo inferior a

dez vezes o número de chamadas feitas pelo

fmincon

nesse problema (~150). Essa

versão modificada do CRSI tem provado excelente taxa de substituição na população,

pois, em média, realiza-se uma iteração para cada menos de dois pontos tentativa

gerados (veja a relação

trial-p

iter

). Os extremos das interpolações quadráticas do

CRSI foram mínimos ou máximos com quase a mesma freqüência, como esperado num

problema geral. Esse aspecto – que não estava muito claro na proposta original de Ali

al.

, 1997 – indica que o CRSI balanceia automaticamente diversificação e

intensificação em seu esquema de busca e acreditamos que essa característica seja a

principal responsável por seu desempenho razoável. As viabilizações forçadas,

empregadas quando alguma coordenada do ponto tentativa é inviável, evita que todas as

suas coordenadas sejam recalculadas, reduzindo-se assim o tempo de CPU do

algoritmo. Esse aspecto é mais importante quando o número de coordenadas do

problema é elevado e sua determinação, custosa. Observe-se também que o ótimo de

algumas variáveis de projeto está sobre um dos seus respectivos limites laterais ou bem

próximo a ele (Tabelas 9 e 10); daí, o esquema de viabilização adotado no CRSI

modificado facilita a convergência da população para essas fronteiras.

A Figura 37 apresenta o histórico de otimização da solução ‘CRSI - 3’ (o mesmo

comportamento geral também é válido para as outras soluções CRSI). Os valores de

função objetivo apresentados correspondem ao melhor ponto da população em cada

iteração (

(

) =

) e incluem a parcela de penalidade. Nas últimas iterações, essa parcela

tende a desaparecer e os valores objetivos tendem para

−

. Os segmentos horizontais

ocorrem quando os pontos tentativa são piores que o atual melhor ponto da população,

ou seja,

– e logo não há melhora no melhor valor objetivo. As mudanças

principais se desenrolam nas primeiras iterações; por outro lado, as últimas iterações

respondem pela contração final da população em torno do ótimo – quando há então

apenas leves modificações no valor da função objetivo.

(a) (b)

Figura 37

−

Histórico de otimização da solução ‘CRSI - 3’ em termos de

(a) iterações e (b) chamadas da função.

Tabela 11 – Comparação do projeto inicial com os melhores resultados de otimização.

Variáveis de Projeto

Grandezas Resultantes

Projeto inicial

(com parametrização)

cubo meio ponta

CRSI - 3

(868

)

cubo meio ponta

fmincon 1

(134

)

cubo meio ponta

(°)

49,3 26,2 17,4 49,5 29,0 18,6 49,4 28,7 18,7

(

−

)

1,61 1,08 0,889 1,65 1,09 0,890 1,61 1,08 0,889

(%)

considerou-se pás retas

3,69 1,77 1,02 3,64 1,67 0,93

(°)

60,0 60,5 60,7

pá

(W) 8761 9178 9218

(W) 215 213 211

(

)

(W) 1443 1011 1011

(%)

80,72 84,71 84,77

(m) 3,70 3,70 3,70

Por fim, a Tabela 11 compara o projeto inicial com os melhores projetos

otimizados produzidos usando-se o SQP-

fmincon

e o CRSI modificado. Como se pode

observar, há uma concordância muito boa entre as soluções ‘CRSI - 3’ e ‘

fmincon

1’.

Esta última tem rendimento ligeiramente superior por causa da característica local da

busca dos métodos de gradiente; de fato, os algoritmos globais, como o CRSA, não tem

0 50 100 150 200 250 300 350 400 450 500

-0.848

-0.846

-0.844

-0.842

-0.84

-0.838

-0.836

-0.834

Iteração

Função objetivo

0 100 200 300 400 500 600 700 800 900

-0.848

-0.846

-0.844

-0.842

-0.84

-0.838

-0.836

-0.834

Número de chamadas do solver

Função objetivo

esse poder de refinamento. De qualquer forma, está bastante aparente que o ótimo

global do problema da turbina foi provavelmente encontrado. Repita-se, porém, que o

SQP-

fmincon

exige um projeto de partida

e que diversos foram testados até que se

tivesse uma certa segurança de que ‘

fmincon

1’ era mesmo a melhor solução. Uma das

principais vantagens do CRSA é exatamente não depender dessas estimativas iniciais,

pois em vez de um único ponto, o CRSA usa uma população inicial.

Como um último comentário, toda essa discussão abre a perspectiva de se usar

em conjunto os métodos local e global (SQP-

fmincon

e CRSI) de pelo menos duas

maneiras. Numa forma seqüencial, inicia-se com o CRSI e, após finalizado o processo

de contração, usa-se a solução encontrada como ponto de partida para o SQP-

fmincon

i.e., fazendo-se

. De fato, isso foi testado algumas vezes para se verificar a

otimalidade dos resultados do CRSI e, em geral, o SQP-

fmincon

acabou produzindo

uma solução bem próxima do ponto de partida, apenas refinando um pouco o resultado

do CRSI. A outra maneira de empregar ambos os métodos é num modo híbrido: toda

vez que o ponto tentativa for o melhor na população (

), “dispara-se” o

fmincon

desse ponto (

) e usa-se a solução local encontrada para substituir o melhor ponto

da população. Embora as buscas se tornassem provavelmente menos diversificadas,

diminuindo-se a chance de encontrar o ótimo global, tal procedimento pode acelerar o

processo de convergência, economizando significativamente o número de chamadas do

solver

do escoamento.

Capítulo 7

CONCLUSÕES E SUGESTÕES

7.1 CONCLUSÕES

Embora se trate de máquinas centenárias, o projeto de turbinas hidráulicas ainda

depende consideravelmente de informações empíricas acumuladas ao longo do tempo.

Nas últimas décadas, porém, as ferramentas de simulação numérica de escoamentos

tornaram-se indispensáveis aos engenheiros, abrindo perspectivas em relação a novas

metodologias de projeto. Neste trabalho, técnicas de otimização são usadas em uma

metodologia de projeto intermediária, de baixo custo computacional, para turbinas

hidráulicas axiais. O nível de sofisticação pretendido enquadra-se na chamada

otimização conceptual. Nesse estágio do dimensionamento da turbina, o objetivo é

determinar os principais parâmetros geométricos do distribuidor e do rotor, tendo-se em

vista a formação de um campo de escoamento favorável ao melhor desempenho da

máquina.

A geometria considerada foi a estritamente necessária para o cálculo dos perfis

de velocidade nas seções de entrada e de saída do distribuidor e do rotor. Como

variáveis de projeto, escolheram-se o ângulo de saída do distribuidor e os ângulos de

montagem, os arqueamentos e os comprimentos dos perfis do rotor, tratados apenas

como linhas médias em arcos de circunferência. Essa seleção de parâmetros permite ao

programa de otimização alterar a quantidade de movimento angular do fluido que o

distribuidor descarrega para o rotor, buscando-se um valor tal que favoreça as melhores

condições de absorção da energia pelas pás, identificando-se as melhorias de

desempenho na fase inicial de projeto. Usaram-se também parametrizações parabólicas,

baseadas nas estações do cubo, meio e ponta das pás, para facilitar a descrição e as

alterações da geometria durante o processo de otimização, assim como para se reduzir o

número de variáveis de projeto.

O cálculo dos perfis de velocidade na saída do distribuidor e do rotor é realizado

considerando-se corretamente a condição de equilíbrio radial das máquinas de fluxo

axiais. Com isso, o campo de velocidade assim determinado acha-se fisicamente

consistente com os princípios de conservação de massa, energia e quantidade de

movimento.

Os desvios e as perdas hidráulicas são calculados usando-se correlações

empíricas da literatura. Buscou-se contemplar as parcelas mais importantes das perdas

que ocorrem em turbinas hidráulicas axiais. O aspecto fundamental é que a maneira

como cada parâmetro de projeto varia as perdas pôde ser evidenciada pelo modelo,

orientando-se a busca em direção a uma geometria de boas características, ou seja, com

perdas mínimas.

De fato, as comparações dos resultados de otimização com um projeto já

existente de turbina hélice tubular mostram que a modelagem de perdas adotada,

juntamente com a parametrização da geometria e o cálculo do equilíbrio radial,

realmente conduzem a tendências de projeto corretas para turbinas hidráulicas axiais e

esse é o êxito maior da metodologia desenvolvida.

Dois métodos de otimização, um local e outro global, foram usados nas buscas

automáticas por geometrias de projeto. Para o método local

−

baseado em gradiente

−

escolheu-se a Programação Quadrática Seqüencial (SQP) implementada na função

fmincon

do MatLab

. Essa técnica, consagrada na programação não-linear com

restrições, mostrou-se conveniente no problema da turbina, apresentando rápida

convergência e necessitando de relativamente poucas avaliações do

solver

/função

objetivo. No entanto, procura apenas mínimos locais e necessita de uma estimativa

inicial, tornando o êxito do processo de otimização bastante dependente desse ponto de

partida. Essas dificuldades são contornadas pelos Algoritmos de Busca Aleatória

•

Controlada (CRSA), o método de busca global escolhido. Trata-se de um método

populacional no qual, partindo-se de uma população inicial – em vez de um único ponto

inicial

−

, executa-se um processo iterativo de convergência dessa população em direção

ao ótimo global, por meios puramente heurísticos.

Uma versão de CRSA foi proposta a partir de algumas modificações sobre um

algoritmo já existente, o CRSI. Aparentemente, os autores do CRSI original não

reconheceram um aspecto do esquema de busca via interpolações quadráticas, que é a

chance de os extremos dessas parábolas serem máximos em vez de apenas mínimos,

como eles sugeriram. De fato, situações de máximo ocorrem freqüentemente, o que

torna as buscas mais exploratórias – e não perturba a contração da população quando

necessário.

A solução “global” obtida usando-se o CRSI está em muito boa concordância

com a melhor solução local encontrada pelo SQP-

fmincon

, de modo que o projeto

ótimo global foi provavelmente obtido.

Em princípio, o projeto otimizado de outros tipos de turbomáquinas axiais –

bombas, ventiladores, turbinas a gás, turbocompressores – pode ser realizado com base

nesta mesma metodologia desenvolvida para turbinas hidráulicas axiais.

7.2 SUGESTÕES

Incluir efeitos de cavitação na formulação do problema de otimização

. De

fato, esse aspecto é básico no projeto de turbinas hidráulicas, foi ignorado até

aqui, mas tem forte influência sobre a forma final dos perfis das pás,

especialmente na definição de espessuras (não consideradas neste trabalho) e

razões corda-passo (

) – estas últimas foram sempre diminuídas até o limite

inferior durante o processo de otimização justamente por não se considerar

questões de descolamento e cavitação. Poderia ser feita a minimização de um

coeficiente de cavitação, como um objetivo adicional do problema, ou então

imporem-se restrições quanto à pressão mínima na superfície das pás. Neste

caso, seria preciso estabelecer e incorporar ao

solver

um módulo para a

distribuição de pressão ao redor dos perfis – talvez um cálculo potencial

•

fosse satisfatório – e a espessura também seria tratada como variável de

projeto.

Usar a classificação de Oyama et al. (2005)

para tratar eficientemente as

restrições e/ou os vários objetivos originados na sugestão anterior. Isso

evitaria dificuldades inerentes ao método de penalização que foi empregado

neste trabalho.

Testar outras formas de parametrização da geometria

. Certamente há outras

funções mais apropriadas para se construir as parametrizações geométricas

do que as simples funções parabólicas usadas neste trabalho. Por exemplo,

poderiam ser investigadas funções bi-lineares, tri-lineares, exponenciais,

hiperbólicas e combinações de funções, todas com fortes possibilidades de

representar satisfatoriamente a geometria do rotor. Além disso, as

parametrizações podem ser baseadas em outras estações além das três usadas

neste trabalho – cubo, meio e ponta. Estudos desse tipo indicariam as formas

mais apropriadas para se descrever a geometria da turbina, sendo importantes

também em problemas tratados via CFD (Drtina e Sallaberger, 1999; Lipej,

2004; Oyama

et al

., 2005).

Tratar corretamente a questão do choque na entrada do rotor

. Considerou-se

como condição sem choque a de incidência nula, o que seria verdadeiro para

um número infinito de pás. Com

pá

finito, a condição sem choque acha-se

com uma pequena incidência diferente de zero e esse detalhe é importante

para que se possa garantir um campo de velocidade favorável na entrada do

rotor. Lembre-se que uma das principais ações do otimizador foi justamente

buscar a condição sem choque.

Testar outros modelos de perdas e desvios

. Sendo escassas as publicações

referentes a turbinas hidráulicas, usou-se correlações empíricas disponíveis

na literatura de turbinas a gás, tendo-se em vista as restrições quanto à sua

aplicação. Seguindo essa linha, outros sistemas de perdas poderiam ser

testados na modelagem da turbina hidráulica axial, como os de Kacker e

Okapuu (1982), Park e Chung (1992), Yoon

et al

. (1998).

Levantar experimentalmente correlações de perdas e desvios para turbinas

hidráulicas axiais

. Esse tipo de investigação raramente aparece na literatura e

seguramente seria mais apropriado a este trabalho do que adaptar resultados

referentes a turbinas a gás.

•

Partir para solver’s mais sofisticados

. Um nível de projeto mais avançado

requer informações mais detalhadas a respeito da geometria e do campo de

escoamento, informações essas que as simplificações usadas neste trabalho já

não são capazes de fornecer. Cálculos potenciais com correção de efeito de

camada-limite ou mesmo programas do tipo 3D-Euler seriam próximas

alternativas. Inevitavelmente, o custo computacional cresce sensivelmente

com essa sofisticação do

solver

. No trabalho apresentado, a meta era uma

metodologia de projeto intermediária, de baixo custo computacional, viável

de se perfazer rapidamente em um único computador (PC) e que precedesse a

análises mais sofisticadas para já adiantar resultados básicos importantes.

Determinar detalhadamente os perfis das pás usando-se uma técnica de

projeto inverso otimizado

. Nessa etapa mais avançada, complementar à

metodologia apresentada, as distribuições de velocidade inicialmente

determinadas na otimização conceptual são preservadas ou levemente

alteradas ao mesmo tempo em que os perfis das pás são projetados de modo a

garantir também mínima ocorrência de cavitação.

Aprimorar o Algoritmo de Busca Aleatória Controlada

. Ficou claro que os

CRSA’s ainda estão abertos a modificações para melhorar seu desempenho,

i.e., sua velocidade de convergência e sua robustez em encontrar mínimos

globais. Durante a execução deste trabalho, estudos paralelos levaram a

versões mais promissoras que o CRSI proposto. De fato, os CRSA’s são

menos difundidos nas aplicações de engenharia.

Testar outras técnicas de otimização global no problema da turbina

, por

exemplo, algoritmos genéticos e de evolução diferencial, mais populares na

comunidade de engenharia. Seus resultados poderiam ser confrontados com

os do CRSA analisando-se históricos de otimização, como os da Figura 37.

Construir a turbina do projeto otimizado, testá-la e compará-la com o

projeto inicial

. Essa seria a prova de fogo! Veríamos se a metodologia de

projeto desenvolvida leva realmente a parâmetros de desempenho ótimos.

Esse estudo seria provavelmente bastante útil para se apontar as deficiências

da modelagem exposta em predizer o comportamento das perdas.

REFERÊNCIAS BIBLIOGRÁFICAS

ALI, M. M., TÖRN, A. e VIITANEN, S. (1997), “A Numerical Comparison of Some

Modified Controlled Random Search Algorithms”, Journal of Global Optimization, Vol.

11, pp. 377-385.

ALI, M. M. e TÖRN, A. (2004), “Population set-based global optimization algorithms: some

modifications and numerical studies” Computers & Operations Research, Vol. 31, pp.

1703-1725.

ALONSO, L. A. P., EL-JAMAL, R. M. e TEIXEIRA, J. C. (1990), “Álgebra II”, Editora

Anglo.

BRAN, R. e SOUZA, Z. (1979), “Máquinas de Fluxo”, Livros Técnicos e Científicos

Editora.

CORDIER, O. (1955), “Ähnlichkeitsbetrachtung bei Strömungsmaschinen”, VDI-Zeitschrift,

Vol. 97, No. 34, pp. 1233-1234.

DANSIE, J. (1996), “Turbines in the 90’s”, International Water Power & Dam Construction,

August, pp. 26-28.

DENTON, J. D. (1993), “Loss Mechanisms in Turbomachines”, ASME Journal of

Turbomachinery, Vol. 115, pp. 621-656.

“DESIGN BY NUMBERS” (1998), International Water Power & Dam Construction, March,

pp. 20-22.

DIXON, L. C. e SZEGÖ, G. P. (1978), “Towards Global Optimisation 2”, North-Holland

Publishing Company.

DRTINA, P. e SALLABERGER, M. (1999), “Hydraulic turbines – basic principles and

state-of-the-art computational fluid dynamics applications”, Proceedings of the

Institution of Mechanical Engineers, Vol. 213, part C, pp. 85-102.

FOX, R. W. e MCDONALD, A. T. (1998), “Introdução à Mecânica dos Fluidos”, Livros

Técnicos e Científicos Editora.

HINDLEY, M. (1996), “Buoyant year for turbines”, International Water Power & Dam

Construction, February, pp. 29.

HORLOCK, J. H. (1973), “Axial Flow Turbines”, Robert E. Krieger Publishing Company.

KACKER, S. C. e OKAPUU, U. (1982), “A Mean Line Prediction Method for Axial Flow

Turbine Efficiency”, ASME Journal of Engineering for Power, Vol. 104, pp. 111-119.

KARAMCHETI, K. (1980), “Principles of Ideal-Fluid Aerodynamics”, Robert E. Krieger

Plublishing Company.

LEE, C. e CHUNG, M. K. (1991), “Secondary Flow Loss and Deviation Models for

Through-Flow Analysis of Axial-Flow Turbomachinery”, Mechanics Research

Communications, Vol. 18(6), pp. 403-408.

LIPEJ, A. (2004), “Optimization method for the design of axial hydraulic turbines”,

Proceedings of the Institution of Mechanical Engineers, Vol. 218, part A, pp. 43-50.

MANZANARES FILHO, N., MOINO, C. A. A. e JORGE, A. B. (2005), “An Improved

Controlled Random Search Algorithm for Inverse Airfoil Cascade Design”, WCSMO,

World Congresses of Structural and Multidisciplinary Optimization, paper No. 4451.

MOURA, M. D. de e BRASIL JÚNIOR, A. C. P. (2005), “Assessment of Turbulent

Modelling for CFD Simulation in Hydroturbines: Draft Tubes”, COBEM, 18

International Congress of Mechanical Engineering, ABCM, Ouro Preto – MG, paper

No. 1920.

NASH, S. G. e SOFER, A. (1996), “Linear and Nonlinear Programming”, McGraw-Hill.

OH, H. W. e KIM, K.-Y. (2001), “Conceptual design optimization of mixed-flow pump

impellers using mean streamline analysis”, Proceedings of the Institution of Mechanical

Engineers, Vol. 215, part A, pp. 133-138.

OSTERWALDER, J. (1960), “Flow measurements on models as a means for determining

the loss distribution in Kaplan turbines”, Escher Wyss News, No. 1, 2 e 3.

OYAMA, A., FUJII, K., SHIMOYAMA, K. e LIOU, M.-S. (2005), “Pareto-Optimality-

Based Constraint-Handling Technique and Its Application to Compressor Design”, 17

AIAA CFD Conference, paper AIAA2005-4983.

PARK, H. D. e CHUNG, M. K. (1992), “Refinement of Spanwise Distribution Models of

Deviation Angle and Secondary Loss for Axial Flow Turbine”, Mechanics Research

Communications, Vol. 19(5), pp. 449-455.

PENG, G., CAO, S., ISHIZUKA, M. e HAYAMA, S. (2002a), “Design optimization of

axial flow hydraulic turbine runner: Part I – an improved Q3D inverse method”,

International Journal for Numerical Methods in Fluids, Vol. 39, pp. 517-531.

PENG, G., CAO, S., ISHIZUKA, M. e HAYAMA, S. (2002b), “Design optimization of

axial flow hydraulic turbine runner: Part II – multi-objective constrained optimization

method”, International Journal for Numerical Methods in Fluids, Vol. 39, pp. 533-548.

PFLEIDERER, C. e PETERMANN, H. (1979), “Máquinas de Fluxo”, Livros Técnicos e

Científicos Editora S.A.

PORTO, M. A. A., GUIMARÃES, A. P. B., ROGAR, M. M., ALMEIDA, J. R. C.,

CASTANHO JÚNIOR, C., CAVALCANTI, M. C. R., FONSECA, D. C. e PINTO

JÚNIOR, J. B. (2005), “The Madeira hydroelectric complex – regional integration and

environmental sustainability using bulb type turbines”, PCH Notícias & SHP News, No.

27, pp. 8-12.

PRICE, W. L. (1977), “A controlled random search procedure for global optimisation”,

Computer Journal, Vol. 20(4), pp. 367-370.

QUANTZ, L. (1976), “Motores hidráulicos”, Editorial Gustavo Gili, S.A.

RAABE, J. (1985), “Hydro Power”, VDI Verlagembh.

SCHWEIGER, F. e GREGORI, J. (1988), “Developments in the design of bulb turbines”,

Water Power & Dam Construction, September, pp. 12-15.

SCHWEIGER, F. e GREGORI, J. (1989), “Design of Large Hydraulic Turbines”, IAHR,

International Symposium on Large Hydraulic Machinery & Associated Equipments,

Beijing, pp. 155-164.

SILVA, E. F. e MIYAZIMA, A. T. (2000), “Cálculo Numérico – Sua Implementação via

Microcomputadores”, apostila de uso interno da UNIFEI, Itajubá, MG.

SITE DA ANEEL: http:\\www.aneel.gov.br – Acessado em Março de 2006.

SITE DO CERPCH: http:\\www.cerpch.unifei.edu.br – Acessado em Março de 2006.

SOUZA JÚNIOR, F., SILVA, R. J., MANZANARES FILHO, N., BARBOSA, J. R. e

TOMITA, J. T. (2005), “Design Point Efficiency Optimization of a Multi-Stage Axial-

Flow Compressor for Aero Application applying a Specially Developed Computer

Code”, COBEM, 18

International Congress of Mechanical Engineering, ABCM, Ouro

Preto – MG, paper No. 1802.

SOUZA, Z. (1989), “Relatórios de teste e operação da turbina hidráulica tipo tubo da MEP”,

Relatório Interno, UNIFEI, Itajubá, MG.

SUBRAHMANYAM, K. S. (1982), “Power Development in India and the Role of Low

Head/Small Scale Hydro”, Fuji Electric Review, Vol. 28, No. 3, pp. 87-88.

SULLEREY, R. K. e KUMAR, S. (1984), “A Study of Axial Turbine Loss Models in a

Streamline Curvature Computing Scheme”, ASME Journal of Engineering for Gas

Turbine and Power, Vol. 106, pp. 591-597.

UEDA, T. (1982), “Improvement of Hydraulic Turbine Efficiency”, Fuji Electric Review,

Vol. 28, No. 2, pp. 34-40.

VIVIER, L. (1966), “Turbines Hydrauliques et leur régulation”, Éditions Albin Michel.

YOON, E. S., KIM, B. N. e CHUNG, M. K. (1998), “Modeling of three dimensional

unsteady flow effects in axial flow turbine rotors”, Mechanics Research

Communications, Vol. 25, No. 1, pp. 15-24.

WILSON, R. B. (1963), “A Simplicial Algorithm for Concave Programming”, Ph.D. Thesis,

Harvard University, Cambridge, MA.

Livros Grátis
( http://www.livrosgratis.com.br )
 
Milhares de Livros para Download:
 
Baixar livros de Administração
Baixar livros de Agronomia
Baixar livros de Arquitetura
Baixar livros de Artes
Baixar livros de Astronomia
Baixar livros de Biologia Geral
Baixar livros de Ciência da Computação
Baixar livros de Ciência da Informação
Baixar livros de Ciência Política
Baixar livros de Ciências da Saúde
Baixar livros de Comunicação
Baixar livros do Conselho Nacional de Educação - CNE
Baixar livros de Defesa civil
Baixar livros de Direito
Baixar livros de Direitos humanos
Baixar livros de Economia
Baixar livros de Economia Doméstica
Baixar livros de Educação
Baixar livros de Educação - Trânsito
Baixar livros de Educação Física
Baixar livros de Engenharia Aeroespacial
Baixar livros de Farmácia
Baixar livros de Filosofia
Baixar livros de Física
Baixar livros de Geociências
Baixar livros de Geografia
Baixar livros de História
Baixar livros de Línguas

Baixar livros de Literatura
Baixar livros de Literatura de Cordel
Baixar livros de Literatura Infantil
Baixar livros de Matemática
Baixar livros de Medicina
Baixar livros de Medicina Veterinária
Baixar livros de Meio Ambiente
Baixar livros de Meteorologia
Baixar Monografias e TCC
Baixar livros Multidisciplinar
Baixar livros de Música
Baixar livros de Psicologia
Baixar livros de Química
Baixar livros de Saúde Coletiva
Baixar livros de Serviço Social
Baixar livros de Sociologia
Baixar livros de Teologia
Baixar livros de Trabalho
Baixar livros de Turismo
 
 