UNIVERSIDADE ESTADUAL PAULISTA “JÚLIO DE MESQUITA FILHO” FACULDADE DE ENGENHARIA CÂMPUS DE ILHA SOLTEIRA RAFAEL VANIN INSERÇÃO DE DISPOSITIVOS PASSIVOS PARA A MITIGAÇÃO DO VÓRTICE DE CORDA EM UMA TURBINA HIDRÁULICA DO TIPO FRANCIS Ilha Solteira 2022 PROGRAMA DE PÓS-GRADUAÇÃO EM ENGENHARIA MECÂNICA RAFAEL VANIN INSERÇÃO DE DISPOSITIVOS PASSIVOS PARA A MITIGAÇÃO DO VÓRTICE DE CORDA EM UMA TURBINA HIDRÁULICA DO TIPO FRANCIS Dissertação apresentada à Faculdade de Engenharia de Ilha Solteira – UNESP como parte dos requisitos para obtenção do título de Mestre em Engenharia Mecânica. Especialidade: Ciências Térmicas. Nome do orientador Prof. Dr. Leandro Oliveira Salviano Nome do coorientador Prof. Dr. Aluísio Viais Pantaleão Ilha Solteira 2022 DEDICATÓRIA Dedico este projeto de pesquisa à minha família que sempre me apoia para que eu possa alcançar meus objetivos. AGRADECIMENTOS Agradeço primeiramente a Deus por permitir a oportunidade de desenvolver esse projeto com dedicação. Agradeço aos meus pais, César e Mariazinha, e meu irmão, Michel, pelo apoio incondicional para realização deste sonho. Agradeço ao meu orientador, Leandro, e coorientador, Aluísio, pelos grandes ensinamentos, direcionamento e todo apoio no projeto de pesquisa. Agradeço à banca de defesa pela participação e contribuição para melhoria do projeto. Agradeço aos demais professores do departamento pelo aprendizado que tive com suas respectivas disciplinas. Agradeço aos meus colegas de pesquisa, Daniel e Henrique, que contribuíram muito para o desenvolvimento do projeto, além da parceria e amizade fora de sala de aula. Agradeço ao meu amigo Arthur M.V. pela motivação ao estudo do mestrado. Agradeço à todos os amigos de mestrado e doutorado que me acolheram e me tornaram parte da equipe nas horas de estudos, me incentivando e me ensinando sempre, e também pela parceria e boas lembranças que levarei comigo. O presente trabalho foi realizado com apoio da Coordenação de Aperfeiçoamento de Pessoal de Nível Superior - Brasil (CAPES) - Código de Financiamento 001. “Se eu vi mais longe, foi por estar sobre ombros de gigantes.”. Isaac Newton. (Trecho de uma carta de Isaac Newton para Robert Hooke, 5 de Fevereiro de 1676) RESUMO A maior parte da energia elétrica gerada no Brasil é proveniente das hidrelétricas, reconhecida por ser uma fonte de energia renovável e não intermitente que pode adequar o despacho conforme a demanda da rede elétrica. A demanda por fontes de energias renováveis vem substancialmente crescendo nos últimos anos em decorrência, entre outros, de aspectos econômicos e ambientais das fontes convencionais, como o petróleo e seus derivados. Porém, algumas destas energias renováveis, como a energia solar e eólica, são consideradas como fontes de energias intermitentes, podendo gerar instabilidades na rede elétrica de distribuição. Desta forma, as hidrelétricas podem operar modulando seu ponto de operação para acomodar estas variações, evitando a falta ou a sobrecarga no fornecimento de eletricidade. Neste sentido, o projeto de turbinas hidráulicas através da modelagem computacional permite a avaliação da dinâmica do escoamento sob várias condições de operação, identificando, ainda na fase de projeto, fenômenos como o vórtice de corda no tubo de descarga, que podem gerar vibrações na turbina, limitando sua faixa de operação e reduzindo sua vida útil. Entretanto, nas análises de pré-projeto em meados da década de 70, as ferramentas de modelagem computacional eram pouco difundidas, ou limitadas, no desenvolvimento das turbinas instaladas nas usinas hidrelétricas do Brasil, quando o desenvolvimento a partir de uma abordagem empírica e/ou analítica era predominante. Desta forma, o objetivo deste trabalho foi investigar a possibilidade de mitigação do vórtice de corda em uma turbina hidráulica benchmarking do tipo Francis, apresentada no Workshop Francis-99, considerando a inserção de dispositivos passivos no tubo de descarga e sobre o perfil da pá do rotor. Para redução no tempo de simulação numérica foi utilizado à técnica de periodicidade do rotor e palhetas direcionadoras móveis. O escoamento é admitido ser incompressível e turbulento. Uma abordagem inicial considerando um regime estacionário é adotado para a definição da geometria do dispositivo passivo mais promissor para mitigar o vórtice de corda, seguido de uma análise em regime transitório. As condições de operação PL (Partial Load), BEP (Best Efficiency Point) e HL (High-Load) foram investigadas. As análises mostraram que dispositivos passivos inseridos sobre às pás do rotor são mais promissores para mitigar o vórtice de corda do que aqueles inseridos no tubo de descarga. Além disso, que ambos dispositivos passivos instalados sobre as pás do rotor ou no tubo de descarga têm significativo impacto no torque e eficiência da turbina hidráulica. Embora a condição PL seja a mais relevante em relação à formação do vórtice de corda para a turbina hidráulica Francis Original (Workshop Francis-99), a inserção de dispositivos passivos também impactou significativamente nas outras condições de operação BEP e HL, indicando, desta forma, que o estudo de dispositivos passivos deve ser realizado mandatoriamente para as três condições operacionais. Contudo, os resultados deste trabalho indica que a inserção de dispositivos passivos é uma técnica promissora para a mitigação do vórtice de corda em turbinas hidráulicas do tipo Francis já instaladas em usinas hidrelétricas, contribuindo, assim, para a revitalização das atuais turbinas instaladas em contrapartida à sua completa substituição. Palavras-chave: turbina hidráulica Francis; vórtice de corda; dispositivos passivos; modelagem computacional. ABSTRACT A large amount of the electricity generated in Brazil comes from hydroelectric plants, recognized for being a renewable and non-intermittent energy source that can adapt the dispatch according to the demand of the electricity grid. The demand for renewable energy sources has been growing substantially in recent years as a result, among others, of economic and environmental aspects of conventional sources, such as oil and its derivatives. However, for some of these renewable energies, such as solar and wind energy, are considered as intermittent energy sources, which can generate instabilities in the distribution grid. In this way, hydroelectric plants can operate by modulating their operating point to accommodate these variations, avoiding lack or overloading the electricity supply. In this sense, the design of hydraulic turbines through computational modeling, allows the evaluation of the flow under various operating conditions, identifying, even in the design phase, phenomena such as the vortex rope in the draft tube, which may generate vibrations in the turbine, limiting its operating range and reducing its useful life. However, in pre- project analyzes in the mid-1970s, computational modeling tools were not widespread, or limited, for the development of turbines installed in hydroelectric plants in Brazil, since the development was done from empirical and/or analytical approaches. Thus, the main goal of the present work is to investigate the possibility of mitigating the vortex rope for a Francis benchmarking hydraulic turbine, presented by Workshop Francis-99, considering the insertion of passive devices at the draft tube and on the rotor blade. To reduce the numerical simulation time-consumption, the rotor and guide-vanes were modeled through periodic technique. The flow is assumed to be incompressible and turbulent. An initial analysis considering a steady- state flow is adopted to define a promising geometry of the passive device to mitigate the vortex rope, which is followed by a transient simulation. The operating conditions PL (Partial Load), BEP (Best Efficiency Point) and HL (High-Load) were investigated. Analyzes have been shown that passive devices inserted on the rotor blades are more promising to mitigate the vortex rope than those inserted in the draft tube. Furthermore, passive devices installed either on the rotor blades or draft tube has a significant impact on the torque and efficiency of the hydraulic turbine. Although the PL condition is the most relevant operating condition related to the formation of the vortex rope for the Original hydraulic turbine (Workshop Francis-99), the insertion of passive devices significantly impacted those other operating conditions BEP and HL, reveling that the study of passive devices must be carry-out mandatory for the three operational conditions. Overall, the results of the present work indicated that the insertion of passive devices is a promising technique for the mitigation of the vortex rope for those Francis hydraulic turbines under operation in current hydroelectric plants, contributing to the revitalization of the hydraulic turbines already installed as a counterpart to its complete replacement. Keywords: Francis hydraulic turbine; vortex rope; passive devices; computational modeling. LISTA DE FIGURAS Figura 1: Geração de eletricidade a partir de diversas fontes no Brasil. ..................... 1 Figura 2: Esquema de usina hidrelétrica. .................................................................... 2 Figura 3 - Faixa de aplicação da turbina hidráulica do tipo Francis. ............................ 3 Figura 4 - Vista aérea da usina hidrelétrica de Ilha Solteira - SP. ............................... 4 Figura 5: Componentes e esquema de montagem das turbinas Francis. ................... 8 Figura 6: Esquema de turbina Francis em corte. ......................................................... 9 Figura 7: Triângulo de velocidades na saída do rotor de cada regime de operação. 10 Figura 8: Fluxo reverso no tubo de descarga. ........................................................... 11 Figura 9: Esquema de vórtice de corda. .................................................................... 11 Figura 10: Exemplo de vórtice de corda. ................................................................... 12 Figura 11: Método por injeção de ar na turbina. ........................................................ 16 Figura 12: Método por injeção de água na turbina. ................................................... 17 Figura 13: Dispositivo passivo do tipo protuberância no tubo de descarga. .............. 18 Figura 14: Variações de protuberância no tubo de descarga. ................................... 18 Figura 15: Dispositivo passivo por alterações no cone do rotor. ............................... 19 Figura 16: Dispositivo passivo do tipo J-grooves. ..................................................... 20 Figura 17: Dispositivo passivo do tipo FRUCE em vermelho. ................................... 21 Figura 18: Domínio computacional. ........................................................................... 26 Figura 19: Domínio computacional e Interfaces. ....................................................... 26 Figura 20: Domínio computacional da voluta. ........................................................... 27 Figura 21: Domínio computacional das palhetas direcionadoras móveis. ................. 28 Figura 22: Domínio computacional do conjunto de pás do rotor. .............................. 29 Figura 23: Domínio computacional do tubo de descarga. ......................................... 30 Figura 24: Esquema apresentando a definição de ortogonalidade. .......................... 32 Figura 25: Refino de malha do canal do rotor. .......................................................... 33 Figura 26: Refino de malha do tubo de descarga. ..................................................... 35 Figura 27: Malha padrão na condição de operação PL com Q-criterion 300 s-2. ....... 37 Figura 28: Comparação Estacionário e Transitório no Regime PL. ........................... 38 Figura 29: Número de Courant nas interfaces dos domínios computacionais. ......... 39 Figura 30: Critério de convergência do erro residual médio (RMS) da pressão e velocidade. ................................................................................................................ 41 Figura 31: Critério de convergência do erro residual médio (RMS) de SST k-ω. ...... 41 Figura 32: Perfil de convergência da eficiência. ........................................................ 42 Figura 33: Perfil de convergência do torque. ............................................................. 42 Figura 34: Critério de Convergência por análise qualitativa. ..................................... 43 Figura 35: Comparação de resultados numéricos e experimental para o Torque (T). .................................................................................................................................. 45 Figura 36: Comparação de resultados numéricos e experimental para a Eficiência. 46 Figura 37: Dispositivos Passivos no tubo de descarga (conjunto 1). ........................ 48 Figura 38: Dispositivos Passivos no tubo de descarga (conjunto 2). ........................ 49 Figura 39: Q-criterion de 300 s-2 dos dispositivos passivos instalados no tubo de descarga (conjunto 1). ............................................................................................... 51 Figura 40: Q-criterion de 300 s-2 dos dispositivos passivos instalados no tubo de descarga (conjunto 2). ............................................................................................... 52 Figura 41: Dispositivos Passivos (DP) Tipo 1, do tipo barreira no rotor. ................... 53 Figura 42: Perfis dos dispositivos passivos do Tipo 1 sobre o rotor. ......................... 54 Figura 43: Q-criterion de 300 s-2 dos DP’s do tipo barreira no rotor. ......................... 55 Figura 44: Detalhe Dispositivo Passivo Tipo 2, do tipo barreira no rotor. .................. 56 Figura 45: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2. ............................... 57 Figura 46: Detalhe Dispositivo Passivo Tipo 2, do tipo barreira no rotor com parâmetros. ............................................................................................................... 58 Figura 47: Comparação do perfil de velocidade com e sem barreira. ....................... 59 Figura 48: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=25 mm, na condição de operação PL. ......................................................................................... 60 Figura 49: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=12,5 mm, na condição de operação PL. ......................................................................................... 61 Figura 50: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=25 mm, na condição de operação BEP. ...................................................................................... 62 Figura 51: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=12,5 mm, na condição de operação BEP. ...................................................................................... 63 Figura 52: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=25 mm, na condição de operação HL. ........................................................................................ 64 Figura 53: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=12,5 mm, na condição de operação HL. ........................................................................................ 65 Figura 54: Resultado da última volta no regime transitório do RU_25_60, na condição de operação PL. ......................................................................................... 68 Figura 55: Resultado da última volta no regime transitório do RU_25_60, na condição de operação BEP. ...................................................................................... 69 Figura 56: Resultado da última volta no regime transitório do RU_25_60, na condição de operação HL. ........................................................................................ 70 LISTA DE TABELAS Tabela 1: Parâmetros de operação da turbina hidráulica Francis a partir do segundo Workshop Francis-99 [9]. .......................................................................................... 24 Tabela 2: Condições de contorno no domínio computacional da voluta. .................. 27 Tabela 3: Condições de contorno no domínio computacional das palhetas direcionadoras móveis. ............................................................................................. 28 Tabela 4: Condições de contorno no domínio computacional do conjunto de pás do rotor. .......................................................................................................................... 29 Tabela 5: Condições de contorno no domínio computacional do tubo de descarga. 30 Tabela 6: Características da discretização do rotor. ................................................. 33 Tabela 7: Propriedades do refino de malha do tubo de descarga. ............................ 34 Tabela 8: Comparação entre hipóteses de regimes estacionário e transitório. ......... 36 Tabela 9: Informações das malhas para o GCI. ........................................................ 40 Tabela 10: Análise de incerteza devido a densidade de malha para o torque e eficiência. .................................................................................................................. 44 Tabela 11: Comparação entre os resultados a partir das malhas N1, N2 e N3 e Workshop Francis-99. ............................................................................................... 44 Tabela 12: Resultados das simulações com DP no tubo de descarga. ..................... 50 Tabela 13: Resultados das simulações com DP do tipo barreira no rotor. ................ 54 Tabela 14: Parâmetros de variação dos dispositivos passivos. ................................ 59 Tabela 15: Resultados de Torque e Eficiência dos dispositivos passivos do tipo 2. . 66 LISTA DE SIGLAS E ABREVIATURAS BEP Ponto Ótimo de Eficiência CFD Dinâmica dos Fluidos Computacional DP Dispositivos Passivos DT Tubo de Descarga EIA Energy Information Administration FRUCE Extensão do cone de rotação livre GCI Índice de Convergência de Malha GGI Interface de Malha Geral HL Carga Alta IEC International Electrotechnical Commission LES Large Eddy Simulation NTNU Universidade Norueguesa de Ciência e Tecnologia NVKS Centro de Energia Hidrelétrica Norueguesa PL Carga Parcial RANS Reynolds-Averaged Navier-Stokes RMS Raiz Quadrática Média RU Rotor SH Capa do rotor SST Shear Stress Transport LISTA DE SÍMBOLOS LATINOS 𝐴𝑖 ⃗⃗ ⃗ Vetor normal da face 𝐶𝑖 ⃗⃗ ⃗ Vetor do centroide da célula ao centroide das células adjacentes 𝑓𝑖⃗⃗ Vetor do centroide da célula à cada centroide das faces da célula 𝑔 Aceleração gravitacional [m/s²] H Altura líquida [m] I Tensor identidade Lc Comprimento do núcleo do jato [mm] Lj Comprimento do jato efetivo [mm] N1, N2, N3 Malha Refinada, Intermediária e Grosseira, respectivamente. P Pressão [Pa] Q Q-criterion [s-2] 𝑄 ̇ Vazão Volumétrica [m³/s] r Distância do fluido até o centro de giro da rotação [m] S Tensor taxa de deformação [s-1] 𝑡 Tempo [s] 𝑇 Torque [N.m] 𝑈,V,W Velocidades do fluido [m/s] 𝑦 Distância do primeiro elemento de malha até a parede [m] z Diferença de altura entre a entrada da voluta e a saída do tubo de descarga [m] LISTA DE SÍMBOLOS GREGOS Ω Tensor taxa de rotação [s-1] 𝜇 Viscosidade dinâmica [kg/m.s] 𝜂 Eficiência da turbina [%] 𝜌 Massa específica [kg/m³] θ Ângulo [°] 𝜔 Velocidade angular [rad/s] LISTA DE SUBÍNDICES abs Absoluto 𝑖, 𝑗, 𝑘 Notação vetorial 1 Entrada 2 Saída SUMÁRIO 1 INTRODUÇÃO ............................................................................................... 1 2 REVISÃO BIBLIOGRÁFICA .......................................................................... 7 2.1 TURBINA HIDRÁULICA DO TIPO FRANCIS ................................................ 7 2.2 VÓRTICE DE CORDA ................................................................................. 10 2.3 IDENTIFICAÇÃO DO VÓRTICE DE CORDA .............................................. 12 2.4 TÉCNICAS PARA MITIGAÇÃO DO VÓRTICE DE CORDA ........................ 15 2.4.1 Método com injeção de fluidos ................................................................. 15 2.4.2 Método por dispositivos passivos ............................................................ 17 3 METODOLOGIA NUMÉRICA ...................................................................... 22 3.1 EQUACIONAMENTO DE EFICIÊNCIA ........................................................ 23 3.2 DOMÍNIO COMPUTACIONAL E CONDIÇÕES DE CONTORNO ............... 25 3.3 DISCRETIZAÇÃO DO MODELO COMPUTACIONAL ................................. 31 3.3.1 Discretização da geometria do Rotor ....................................................... 31 3.3.2 Discretização da geometria do Tubo de Descarga .................................. 34 3.4 ANÁLISE PRELIMINAR ............................................................................... 35 3.5 ESTUDO DA DENSIDADE DE MALHA ....................................................... 40 3.6 VALIDAÇÃO NUMÉRICA ............................................................................. 44 4 RESULTADOS E DISCUSSÕES ................................................................. 47 4.1 RESULTADOS PARA DISPOSITIVOS PASSIVOS INSTALADOS NO TUBO DE DESCARGA ........................................................................................................ 47 4.2 RESULTADOS PARA DISPOSITIVOS PASSIVOS SOBRE O ROTOR TIPO 1 ..................................................................................................................... 53 4.3 RESULTADOS PARA OS DISPOSITIVOS PASSIVOS SOBRE O ROTOR TIPO 2 ..................................................................................................................... 56 4.4 RESULTADOS PARA O DISPOSITIVO PASSIVO TIPO 2 COM VARIAÇÃO NOS PARÂMETROS GEOMÉTRICOS ..................................................................... 58 4.5 RESULTADOS DO RU_25_60 NO REGIME TRANSITÓRIO ...................... 67 5 CONCLUSÕES ............................................................................................ 71 6 RECOMENDAÇÕES PARA TRABALHOS FUTUROS ............................... 73 REFERÊNCIAS ............................................................................................ 74 1 1 INTRODUÇÃO Com o aquecimento global e a necessidade de redução do uso de energia proveniente do petróleo, diversos estudos técnicos têm sido conduzidos sobre a aplicação de energias alternativas ou renováveis para a geração de eletricidade em grande escala. Entre as principais formas de energias, destacam-se as fontes renováveis provenientes de matéria orgânica, seja animal e/ou vegetal (biomassa), força dos ventos (eólica), captação da luz do sol (solar) e hidráulica (37). A Figura 1 apresenta um histórico de geração de eletricidade no Brasil entre os anos de 2000 a 2020, obtidos a partir da agência americana Energy Information Administration (EIA), a qual é responsável por coleta, análise e disseminação de informações sobre energia no mundo (9). Como pode ser observado, a energia elétrica a partir das usinas hidrelétricas teve uma queda na participação relativa comparando com as demais fontes geradoras, embora permaneça como a principal fonte de geração de eletricidade. Figura 1: Geração de eletricidade a partir de diversas fontes no Brasil. Fonte: EIA (2021). As usinas hidrelétricas são responsáveis pela geração da maior parte da energia elétrica fornecida para a rede de distribuição no Brasil. Além do 0% 10% 20% 30% 40% 50% 60% 70% 80% 90% 100% G e ra ç ã o d e E le tr ic id a d e R e la ti v a Ano Hidrelétrica Solar Eólica Biomassa Combustível Fóssil Nuclear 2 fornecimento de energia com controle de tensão e frequência, a forma construtiva das barragens retem a água em reservatórios abertos que proporcionam o controle parcial da vazão dependendo da ocorrências das chuvas e regime de operação da hidrelétrica, contribuindo, desta forma, para o controle de cheias, irrigação, processamento industrial, recreação e navegação através de eclusas. A Figura 2 mostra um desenho esquemático de uma usina hidrelétrica onde a água é conduzida pelo conduto forçado sendo direcionado para as turbinas, e então retorna ao fluxo do rio passando pelo tubo de descarga. Os geradores são acoplados às turbinas e geram eletricidade para a rede de transmissão. Figura 2: Esquema de usina hidrelétrica. Fonte: Imagem adaptada de Voith (2013). Cada projeto de uma usina hidrelétrica torna-se único devido às características do local, necessitando de uma análise detalhada, ponderando fatores para estimar e avaliar a viabilidade do projeto. A escolha do tipo de turbina a ser utilizada na planta geradora é um fator importante, que depende das especificações de altura de queda d’água e vazão. Henn (17) apresenta a faixa de aplicação de cada tipo de turbina dependendo da altura de queda d’água (em metros) e da vazão 3 (em m³/s), conforme apresentado na Figura 3. As turbinas do tipo Francis se destacam com relação às demais devido a ampla faixa de aplicação. Figura 3 - Faixa de aplicação da turbina hidráulica do tipo Francis. Fonte: Henn (2016). O desenvolvimento de turbinas hidráulicas do tipo Francis começou a partir de 1848, quando James B. Francis, Engenheiro Chefe da Locks and Canals Company em Lowell, no estado de Massachusetts dos Estados Unidos, realizou uma série de testes no intuito de melhorar a validação da eficiência de turbinas hidráulicas. Após diversos testes, em 1920 o projeto básico de uma turbina surgiu e ficou conhecida como turbina hidráulica do tipo Francis, homenageando James pela sua contribuição (26). Esta turbina foi modelada para converter a energia de pressão da entrada em energia cinética no rotor, que quando acoplado ao eixo de um gerador elétrico produz energia elétrica (3). A utilização das turbinas do tipo Francis é amplamente difundida, com aplicações nas usinas hidrelétricas de Itaipu em Foz do Iguaçu – PR (18), Ilha 4 Solteira – SP (5) (vista aérea apresentada na Figura 4), Três Gargantas na China (46), e Tokke na Noruega (36). Figura 4 - Vista aérea da usina hidrelétrica de Ilha Solteira - SP. Fonte: CESP (2015). Porém, como a utilização deste tipo de turbina ocorreu a partir de 1920, existem instalações com projetos muito diferentes quando comparados com projetos mais atuais, que foram desenvolvidos a partir de moderna modelagem computacional. Khalil (20) aponta que as pesquisas correlatas à área da Dinâmica dos Fluidos Computacional (CFD) começaram a ser feitas por volta do início da década de 70, considerando um escoamento de fluido em um tubo, com geometria simples, simétrico, sem recirculação e vórtices, tendo a velocidade e a pressão como parâmetros principais avaliados através da solução das Equações de Navier-Stokes. Em anos seguintes, começaram a surgir modelos de turbulência (20). Entretanto, somente a partir de 1986 as simulações numéricas (CFD) começaram a ser utilizadas intensamente em pesquisas acadêmicas e na indústria. Muitas dificuldades eram reportadas em relação às condições de contorno em superfícies irregulares, uma vez que as malhas empregadas eram do tipo ortogonal simples, o que tornava difícil modelar geometrias complexas. Além disso, devido às limitações computacionais da época, a convergência da solução era lenta e, consequentemente, o tempo computacional muito elevado (20). Atualmente, a velocidade de processamento tem aumentado vertiginosamente, permitindo soluções numéricas mais robustas e aplicação em problemas cada vez mais complexos, como no projeto de turbinas hidráulicas, 5 minimizando eventuais características na dinâmica do escoamento que resultaria em baixa eficiência de conversão de energia. Neste sentido, o Centro de Energia Hidrelétrica Norueguesa (NVKS) tem se tornado um centro de referência mundial para desenvolvimento e educação em tecnologia de energia hidrelétrica, sediado na Universidade Norueguesa de Ciência e Tecnologia (NTNU). Neste centro há cooperação entre universidades, instituições de pesquisa e empresas hidrelétricas. O NVKS produziu uma série de Workshops, chamada Francis-99, sobre um modelo de turbina Benchmarking do tipo Francis, em escala 1:5,1 de uma turbina real utilizada na usina hidrelétrica de Tokke (21). A NTNU (36) disponibiliza em seu site artigos científicos relacionados aos estudos desenvolvidas em turbinas hidráulicas, geometrias, malhas computacionais e informações sobre as condições de operação da turbina Francis-99. A partir destes modelos Benchmarking disponibilizados, é possível desenvolver estudos em turbinas hidráulicas, simulando o equipamento em condições de operação fora do ponto de projeto BEP (sigla do inglês Best Efficiency Point), situação operacional que pode ocorrer quando há a inserção, por exemplo, de novas fontes geradoras de energia na rede de distribuição elétrica. A consequência do acoplamento de diversas fontes de energia é a necessidade da usina hidrelétrica modular e operar fora do ponto de operação ótimo, atuando como um sistema pontualmente regulador, evitando instabilidade da rede que acarretaria em possíveis problemas ao consumidor final (24). Quando a turbina opera em carga parcial, fora do ponto nominal de projeto, surge uma estrutura de vórtice no tubo de descarga gerados por velocidades tangenciais na saída do rotor, conhecido como vórtice de corda (24). Este fenômeno impacta no aumento da vibração no equipamento devido à oscilação desta estrutura de vórtices no tubo de descarga, dificultando a operação da turbina e prejudicando a durabilidade do equipamento e seus componentes (16). Uma interessante alternativa para contribuir para a mitigação deste fenômeno é a inserção de dispositivos passivos sobre o rotor ou no tubo de descarga, a partir da investigação de diferentes geometrias, tamanho e local de instalação (24). Dispositivos passivos têm sido aplicados com sucesso em áreas como a aeronáutica (14) e trocadores de calor (40), promovendo a manipulação da dinâmica do escoamento para atingir objetivos de aumento da sustentação, e intensificação da transferência de calor, por exemplo. 6 Tendo em vista as situações atuais das usinas hidrelétricas, foi abordado neste trabalho, através da Dinâmica dos Fluidos Computacional (CFD), o impacto de diversos tipos de dispositivos passivos no vórtice de corda gerado no tubo de descarga de uma turbina hidráulica do tipo Francis, considerando os três pontos operacionais de vazão de água, conforme abordado no Workshop Francis-99: PL (Carga Parcial), BEP (ponto nominal) e HL (Carga alta). A abordagem numérica será realizada a partir da geometria disponibilizada pelo Workshop Francis-99 da NTNU (36), considerando como estratégia de análise uma abordagem estacionária e transitória a partir de uma geometria periódica do rotor e das palhetas direcionadoras. 7 2 REVISÃO BIBLIOGRÁFICA A energia elétrica proveniente das hidrelétricas ainda é a maior fonte de energia elétrica no Brasil. Entretanto, a inserção de fontes intermitentes, como a energia solar e eólica, é um caminho natural diante do aumento da demanda e dos impactos ambientes associados a outras fontes a partir de combustíveis fósseis, ou até mesmo à resistência de expansão de usinas hidrelétricas de grande porte. A inserção destas fontes de energias intermitentes na rede elétrica interligada pode produzir oscilações na rede de distribuição. Como as turbinas hidráulicas das usinas hidrelétricas têm capacidade de operação flexível, sua aplicação como reguladora da geração elétrica tem sido cada vez mais comum, submetendo a unidade geradora a operar fora da condição nominal de projeto (45). Quando a turbina se encontra em condições de baixa vazão (PL), por exemplo, surge o fenômeno do vórtice de corda no tubo de descarga, impactando a condição de operação segura da turbina, especialmente no que se refere à instabilidade devido à vibração provocada pela interação da dinâmica do escoamento e as paredes do tubo de descarga. As próximas seções mostram os principais sistemas de uma turbina hidráulica do tipo Francis, bem como detalhes do vórtice de corda e as investigações que têm sido realizadas para a sua mitigação. 2.1 TURBINA HIDRÁULICA DO TIPO FRANCIS A Figura 5 apresenta um esquema geral da montagem de turbina hidráulica do tipo Francis, onde a água escoa pelo conduto forçado antes de entrar na caixa espiral (voluta). Enquanto escoa, a água entra pelo pré-distribuidor (palhetas direcionadoras fixas) para direcionar e uniformizar o escoamento na turbina. O controle de vazão é feito através das palhetas direcionadoras por meio de servomotores (palhetas direcionadoras móveis), restringindo a vazão conforme demanda. O escoamento alcança o rotor da turbina promovendo a rotação do eixo, que quando acoplado ao eixo de um gerador elétrico produz eletricidade. A água é direcionada para o tubo de sucção (tubo de descarga) e posteriormente é liberada para o curso do rio. Peças menores, como anéis de união e tampão, são utilizadas para interligar e isolar um mecanismo do outro no conjunto da turbina (18). 8 Figura 5: Componentes e esquema de montagem das turbinas Francis. Fonte: Itaipu (2021). As turbinas Francis são equipamentos do tipo escoamento combinado, onde o fluido entra radialmente pela voluta com relação ao eixo, e sai axialmente pelo tubo de descarga, conforme esquema apresentado na Figura 6. É também considerada uma turbina de reação, por converter parcialmente a energia hidráulica disponível (energia de pressão) em energia cinética antes de entrar no rotor. A energia de pressão na entrada do rotor é maior do que na saída, pois é convertida em energia cinética ao longo da passagem do fluido pelo rotor. A rotação do rotor ocorre devido à força de sustentação gerado pelo escoamento de água sobre o perfil das pás da turbina (23). 9 Figura 6: Esquema de turbina Francis em corte. Fonte: Imagem adaptada da Voith (2021). O tubo de descarga tem um papel importante no conjunto da turbina, que é converter a energia cinética em energia de pressão, podendo recuperar cerca de 80% do total da pressão somente na parte do cone, o que representa aproximadamente 10% do comprimento total do tubo de descarga partindo da entrada do cone (28). O perfil de escoamento neste componente pode variar consideravelmente dependendo das condições de operação, principalmente se estiver operando fora do ponto nominal da turbina, podendo gerar vórtice no tubo de descarga. Quando a turbina opera em condições fora do ponto nominal (BEP), surgem componentes de velocidades tangenciais na saída do rotor (24), e dependendo da frequência de oscilação dessas componentes há a formação do vórtice de corda, que pode causar instabilidades na turbina, fadiga nas pás do rotor e até mesmo redução da vida útil dos componentes da turbina (27). A Figura 7 apresenta o triângulo de velocidade na saída do rotor, onde V representa a velocidade absoluta entre as componentes de velocidade axial (Vm) e tangencial (Vt). U é a velocidade periférica na saída do rotor e W é a componente 10 relativa da velocidade absoluta. Na condição BEP, existe velocidade apenas no sentido axial, sem velocidade tangencial. Já na condição PL (sigla do inglês Partial Load), além de componentes axiais de velocidade, há componentes tangenciais, uma vez que o escoamento de água é menor que na condição BEP. Neste caso, o escoamento é direcionado para as periferias do tubo devido a força centrífuga, de modo que gera vórtices no mesmo sentido de rotação do rotor. Na condição de operação HL (sigla do inglês High Load), existem também componentes axiais e tangenciais de velocidade na saída do rotor, e como o escoamento de água é maior do que na configuração BEP, a componente de velocidade tangencial tem sentido contrário da rotação do rotor, produzindo vórtices que giram em sentido oposto ao do rotor (24). Figura 7: Triângulo de velocidades na saída do rotor de cada regime de operação. Fonte: Autoria Própria. 2.2 VÓRTICE DE CORDA Embora Kumar (24) afirme que o vórtice de corda é um fenômeno que ocorre em turbinas hidráulicas devido à existência de velocidade tangencial do escoamento na saída do rotor, ou seja, apenas na condição PL e HL, gerando um vórtice na saída do rotor e se alongando dentro do tubo de descarga. Porém, na condição BEP também apresenta uma estrutura de vórtice centralizada, devido a dinâmica do escoamento. A intensidade da velocidade tangencial altera a velocidade do escoamento no centro do tubo de descarga gerando um fluxo reverso, conforme apresentado na Figura 8 (49). 11 Figura 8: Fluxo reverso no tubo de descarga. Fonte: Zobeiri (2009). Kumar (24) explica que o fluxo reverso gera uma zona de recirculação criando uma região de estagnação com relação à passagem da água pelo tubo de descarga. O vórtice de corda gira em torno desta região conforme apresentado na Figura 9. Figura 9: Esquema de vórtice de corda. Fonte: Autoria Própria. 12 Um exemplo do vórtice de corda foi identificado por Kobro (21) durante experimentos com uma turbina Francis na NTNU sob o regime PL, conforme mostrado na Figura 10. O vórtice de corda é um fenômeno periódico, e sua geometria e frequência de oscilação dependem do projeto da turbina e do condição de operação (49). Figura 10: Exemplo de vórtice de corda. Fonte: Kobro (2010). 2.3 IDENTIFICAÇÃO DO VÓRTICE DE CORDA A existência de um vórtice de corda é detectada a partir da variação da pressão na parede do tubo de descarga, causada pela oscilação de pressão que o vórtice de corda está provocando. Os efeitos prolongados desse fenômeno no tubo de descarga provocam danos em componentes da turbina e estruturas da usina hidrelétrica, uma vez que há uma distribuição assimétrica de pressão gerando forças radiais rotativas que resultam em vibração, reduzindo assim a vida útil devido à fadiga (8). A presença do fenômeno do vórtice de corda aumenta as perdas hidráulicas da turbina devido a redução da eficiência de conversão de energia que deveria ocorrer no tubo de descarga (41). 13 Desse modo, a identificação deste vórtice de corda a partir de uma modelagem computacional torna-se importante para o estudo de seu comportamento ao longo do tubo de descarga. Zhang (48) realizou revisão dos diversos métodos existentes para a identificação de vórtices em simulações numéricas, apresentando as particularidades de cada método, embora o autor não seja taxativo em afirmar qual é o mais adequado para avaliar vórtices de corda em turbinas Francis. Zhang (48) aborda o método da pressão, que caracteriza os vórtices graficamente como uma superfície de mesmo valor de pressão (termo do inglês Iso-surface Pressure), sendo o vórtice identificado por apresentar uma pressão central muito menor que a pressão ao seu redor. Minakov (33) e Gravilov (13) escolheram este método com valores de referência de pressão que melhor se assemelha com o vórtice de corda dentro do tubo de descarga, porém não reportaram os valores usados. Minakov (33) utiliza esse método para realizar uma análise na qual compara o espectro de frequência de oscilação de pressão no tubo de descarga obtidos de simulações numéricas com experimentos realizados em laboratório, e os resultados foram próximos aos experimentais. Gravilov (13) utiliza o mesmo método para caracterização do comportamento do vórtice de corda durante o acionamento da turbina Francis-99, simulando a abertura das palhetas direcionadoras móveis iniciando no regime de operação PL até a abertura e vazão correspondente ao BEP. Durante o estudo numérico, é evidenciado que o vórtice de corda do regime PL se dissipa, permanecendo apenas uma estrutura central, característica da turbina operando em regime BEP. Outro método para a identificação do vórtice de corda conhecido por Q- criterion é abordado por Zhang (48). Este método caracteriza o vórtice através da segunda invariante do tensor gradiente de velocidade (Q), conforme indicado na Equação (1). Jeong (19) desenvolveu o equacionamento para a obtenção do valor de Q: 𝑄 = 1 2 [|Ω|2 − |𝑆|2] > 0 (1) onde: a norma de Ω (Tensor taxa de rotação) é obtido pela Equação (2): |Ω|2 = tr [ ΩΩ 𝑇 ] (2) 14 E Ω é obtido pela Equação (3): Ω⃗⃗⃗⃗ ⃗⃗⃗⃗ = 1 2 [ (∇ �⃗⃗� )𝑇 − ∇ �⃗⃗� ] (3) Já a norma de S (Tensor taxa de deformação) é obtido pela Equação (4): |𝑆|2 = tr [ 𝑆𝑆𝑇 ] (4) E S é obtido pela Equação (5): 𝑆 = 1 2 [ (∇ �⃗⃗� ) 𝑇 + ∇ �⃗⃗� ] (5) Para uma análise qualitativa do vórtice de corda utilizando o Q-criterion, os autores Tran (43), Krane (22) e Fahlbeck (10), sugerem valores de Q maiores que zero para identificar os vórtices. Tran (43) simulou a turbina Francis-99 completa em condição de operação com modelo de turbulência SST k-ω (sigla do inglês Shear Stress Transport) utilizando o valor de Q igual a 300 s-2. Comparando os resultados de Tran (43) com o resultado experimental obtido por Kobro (21), apresentado na Figura 10, Tran (43) concluiu que os resultados foram qualitativamente satisfatórios. Krane (22) simulou a turbina Francis-99 na condição operacional nominal (BEP), considerando diferentes modelos de turbulência e utilizando o Q-criterion com o valor de 200 s-2. Krane (22) concluiu que o modelo de turbulência com melhor concordância com a abordagem experimental foi o LES (sigla do inglês Large Eddy Simulation). Fahlbeck (10) também simulou a turbina Francis-99 com o modelo de turbulência LES em condição de desligamento, partindo do BEP até o completo fechamento das palhetas direcionadoras móveis. Fahlbeck (10) avaliou o comportamento do vórtice de corda através do método Q-criterion com o valor de 200 s-2, para os ângulos de abertura das palhetas direcionadoras móveis de 9,84° e 4,33°, verificando que houve a ruptura do vórtice de corda principal em vários vórtices menores e dispersos. Contudo, é possível observar que o método do Q-criterion tem sido preferencialmente utilizado para a caracterização do vórtice de corda a partir de 15 simulações numéricas, embora o valor do Q-criterion deva ser adequado em função das características do refino da malha e também da condição operacional da turbina. 2.4 TÉCNICAS PARA MITIGAÇÃO DO VÓRTICE DE CORDA Algumas pesquisas têm sido realizadas com o objetivo de mitigar o fenômeno da formação do vórtice de corda, conforme pode ser visto nos trabalhos de Luo (29) e Qian (38). Duas abordagens têm sido mais investigadas na literatura aberta, a injeção de fluidos e a inserção de dispositivos passivos, os quais são discutidos a seguir. 2.4.1 Método com injeção de fluidos A técnica de mitigação do vórtice de corda em turbina hidráulica do tipo Francis através da injeção de fluido consiste em injetar água ou ar no início do tubo de descarga. Este método requer testes e simulações para definir a quantidade necessária de fluido de injeção em cada condição operacional. Por necessitar de sistema externo com a necessidade de energia, esta abordagem é considerada uma técnica ativa de mitigação. Luo (29) simulou uma turbina do tipo Francis na condição operacional PL, injetando diferentes vazões de ar no tubo de descarga, e concluiu que conforme a vazão do ar foi aumentada a oscilação na pressão no tubo de descarga foi gradualmente reduzida. March (30) realizou três métodos de injeção de ar na região central, periférica e distribuída no início do tubo de descarga, conforme esquema mostrado na Figura 11. March realizou testes sob diferentes condições operacionais para várias proporções de fluido/ar, concluindo que independente do ponto de operação e do método, houve redução de eficiência da turbina hidráulica quando se aumentava a proporção de ar. 16 Figura 11: Método por injeção de ar na turbina. Fonte: March (2011). Foroutan (11) realizou um estudo para controlar a formação do vórtice de corda no tubo de descarga com a injeção de água na região axial do tubo por meio de simulações numéricas, conforme mostrado na Figura 12. Foroutan (11) concluiu que os resultados foram promissores, pois eliminou a região de estagnação no centro do tubo de descarga e, consequentemente, o vórtice de corda, reduzindo o ruído e a oscilação na pressão no tubo de descarga. 17 Figura 12: Método por injeção de água na turbina. Fonte: Foroutan (2014). 2.4.2 Método por dispositivos passivos O método por dispositivos passivos consiste na inserção de geometrias no tubo de descarga ou sobre as superfícies das pás do rotor. Estas geometrias não requerem energia para seu funcionamento e, portanto, são consideradas técnicas passivas. Algumas abordagens vêm sendo desenvolvidas para mitigar o efeito do vórtice de corda a partir da inserção de protuberâncias (35; 39), extensão do cone do rotor (15; 38) e ranhuras no tubo de descarga e/ou no rotor (6; 38). Rheingans (39) realizou experimentos inserindo protuberâncias no tubo de descarga, conforme mostrado na Figura 13, com o objetivo de uniformizar o escoamento e reduzir a oscilação na pressão. Rheingans (39) observou que a instalação das protuberâncias próxima à saída do rotor reduziu a amplitude e a magnitude da oscilação do vórtice de corda, reduzindo também a magnitude da vorticidade na saída do tubo de descarga. O pesquisador concluiu que o efeito da oscilação da potência gerada pela turbina hidráulica sempre é acompanhada pela oscilação de pressão na comporta e oscilação do vórtice de corda, e que a inserção de dispositivos passivos no tubo de descarga reduziu as oscilações de potência. 18 Figura 13: Dispositivo passivo do tipo protuberância no tubo de descarga. Fonte: Rheingans (1940). Nishi (35) também realizou experimentos utilizando protuberâncias com perfis triangular e trapezoidal, conforme mostrado na Figura 14. Nishi (35) observou que a frequência do vórtice aumenta conforme aumenta a altura da protuberância, e concluiu que estes tipos de protuberâncias ampliavam a faixa de operação quando instaladas no início do tubo de descarga de uma turbina Francis. Figura 14: Variações de protuberância no tubo de descarga. Fonte: Nishi (1996). Qian (38) simulou três diferentes modificações no cone do rotor de uma turbina do tipo Francis, considerando a extensão do cone, extensão com grooves (traduzindo do inglês, ranhuras), e extensão com topo arredondado, conforme 19 mostrado na Figura 15. Qian (38) concluiu que as três modificações alteraram positivamente o comportamento do vórtice de corda no tubo de descarga, reduzindo a oscilação na pressão na região de descarga. O autor indica também que o modelo com cone estendido com topo arredondado foi o que mais surtiu efeito sobre o vórtice de corda. Chen (6) realizou simulações para a otimização do formato do dispositivo passivo do tipo J-groove inseridos no tubo de descarga, conforme mostrado na Figura 16. Chen (6) concluiu que a inserção deste tipo de dispositivo passivo reduziu a intensidade do vórtice de corda e a oscilação na pressão no tubo de descarga. Figura 15: Dispositivo passivo por alterações no cone do rotor. Fonte: Qian (2012). 20 Figura 16: Dispositivo passivo do tipo J-grooves. Fonte: Chen (2019). Gogstad (15) realizou experimentos com um dispositivo passivo denominado FRUCE (sigla do inglês Free Rotating Cone Extension), o qual é uma extensão do cone conectado ao cubo do rotor, investigando a variação no comprimento do dispositivo, conforme destacado em vermelho na Figura 17. Gogstad (15) notou que sob a condição operacional PL com o dispositivo maior (Figura 17c), a eficiência foi levemente aumentada, enquanto os demais dispositivos reduziram a eficiência da turbina hidráulica. O autor observou ainda redução na eficiência sob a condição operacional BEP, e que na condição operacional HL houve um leve aumento de eficiência com a utilização do maior dispositivo. Gogstad (15) finaliza concluindo que o dispositivo FRUCE é uma técnica efetiva para manipular as oscilações de pressão no tubo de descarga associados às modificações no vórtice de corda. 21 Figura 17: Dispositivo passivo do tipo FRUCE em vermelho. Fonte: Gogstad (2017). Conforme verificado, as linhas de pesquisas sobre a mitigação do vórtice de corda não são recentes. Rheingans (39) indica que em 1940 percebeu oscilações na potência de uma turbina hidráulica e que através de dispositivos passivos a oscilação do vórtice de corda poderia ser reduzida. A solução mais adequada para a mitigação do vórtice de corda ainda não é conhecida, especialmente devido às características do projeto da turbina e da condição operacional, destacando a necessidade de investigação da mais promissora técnica passiva considerando a geometria, tamanho, quantidade e local de instalação. Além disso, os métodos de dispositivos passivos se tornam ainda mais atrativos por não depender da utilização de energia externa durante o seu uso na turbina. Assim, é evidente que a utilização de uma abordagem numérica para o desenvolvimento de dispositivos passivos aplicados a turbinas hidráulicas do tipo Francis que estejam em operação é relevante para tornar a integração entre uma usina hidrelétrica e fontes de energias intermitentes ainda mais atrativos. 22 3 METODOLOGIA NUMÉRICA A modelagem da dinâmica do escoamento em uma turbina hidráulica do tipo Francis, sob diferente condições operacionais, é desenvolvida através da solução das equações de Navier-Stokes, conforme abordagem utilizada por Minakov (33). Para a solução das equações governantes, o software Ansys CFX 19.0 foi empregado, onde uma abordagem acoplada pressão-velocidade é definida, resolvendo por método iterativo em um único sistema de equações algébricas [38]. Funções de interpolação de segunda ordem foram usadas e a discretização das equações é feita a partir do Método dos Volumes Finitos (34). Foram consideradas as hipóteses de escoamento tridimensional, incompressível, turbulento e em regime transitório, em concordância com Choi (7). Dessa maneira a equação da continuidade e quantidade de movimento podem ser apresentadas pelas Equações (6) e (7), respectivamente, como: ∇ ∙ (ρ �⃗⃗� ) = 0 (6) 𝜕(𝜌 �⃗⃗� ) 𝜕𝑡 + ∇ ∙ (𝜌 �⃗⃗� ⊗ �⃗⃗� ) = −∇ 𝑝 + ∇ ∙ [ 𝜇 (∇ �⃗⃗� + (∇ �⃗⃗� ) 𝑇 − 2 3 𝐼 ∇ ∙ �⃗⃗� )] − 2 𝜌 �⃗⃗� × �⃗⃗� − 𝜌�⃗⃗� × (�⃗⃗� × 𝑟 ) (7) onde: 𝜌 é a massa específica do fluido, em kg/m³, �⃗⃗� é a velocidade de escoamento local do fluido, em m/s, 𝑝 é a pressão no fluido, em Pa, 𝜇 é a viscosidade dinâmica do fluido, em Pa.s, 𝐼 é o tensor identidade, adimensional, �⃗⃗� é a velocidade angular do fluido, em rad/s, 𝑟 é a distância do fluido até o centro de giro da rotação, em m. A caracterização do vórtice de corda indica que uma abordagem considerando um regime transitório é a mais adequada. Entretanto, para um passo- angular de 1° por time-step (traduzindo do inglês, passo de tempo), conforme sugerido por Tran (43), demanda um enorme dispêndio de recurso computacional. A necessidade de uma malha refinada, associado ao passo-angular pode conduzir a um tempo de análise proibitivo considerando os recursos computacionais disponíveis para este trabalho. Além disso, conforme discutido na seção 2.4.2, a definição de um dispositivo passivo para a mitigação do vórtice de corda é uma área em evolução, necessitando muitas análises para a identificação de uma técnica 23 passiva promissora. Desta forma, a estratégia adotada neste trabalho inicialmente considerou a investigação de diversas alternativas de dispositivos passivos sob a hipótese de regime estacionário, observando, naturalmente, as características qualitativas dos vórtices identificados como aquelas mostradas na Figura 10. Assim, o dispositivo passivo inicialmente identificado como promissor para a mitigação do vórtice de corda foi submetido a uma análise transitória para a quantificação de seu efeito sobre a eficiência e torque da turbina hidráulica Francis, bem como sobre a formação do vórtice de corda. Para a modelagem computacional URANS (Unsteady. Reynolds Averaged Navier-Stokes), o modelo de turbulência SST k-ω apresenta características adequadas ao tipo de escoamento em turbinas hidráulicas. Choi (7) realizou a simulação numérica em regime transitório em uma turbina do tipo Francis, sob diferentes pontos operacionais, com o modelo de turbulência SST k-ω, justificando ser um modelo adequado quando o escoamento está sujeito a um gradiente adverso de pressão. Yaping (47) e Laouari (25) também realizaram simulações em turbinas hidráulicas do tipo Francis com o modelo de turbulência SST k-ω. Além destas observações, o modelo SST k-ω também é indicado para identificar regiões suscetíveis ao descolamento da camada limite e grande distorção no escoamento (31). 3.1 EQUACIONAMENTO DE EFICIÊNCIA A partir dos dados divulgados no Primeiro e Segundo Workshop Francis-99 (36), informações importantes para a modelagem computacional e cálculo de eficiência podem ser adquiridos, tais como: pressão, torque, vazão, rotação do rotor, densidade da água, ângulo de abertura das palhetas direcionadoras móveis e altura líquida, conforme norma IEC 60193:1999 (IEC, abreviação do inglês International Electrotechnical Commission) (36). Entretanto, o equacionamento para o cálculo de eficiência e altura líquida disponível do segundo Workshop Francis-99 difere daquele reportado no primeiro Workshop Francis-99. Como o presente trabalho utiliza o segundo Workshop Francis-99 como referência para a modelagem computacional, adotou-se a metodologia descrita neste relatório. A Tabela 1 apresenta os dados para os experimentos no segundo Workshop (36), os quais são utilizados para a 24 realização das simulações numéricas. Aggidis (1) realizou procedimento semelhante ao adotado neste trabalho para o cálculo da eficiência da turbina hidráulica do tipo Francis. Tabela 1: Parâmetros de operação da turbina hidráulica Francis a partir do segundo Workshop Francis-99 (36). Parâmetros PL BEP HL Ângulo da Palheta direcionadora móvel (°) 6,72 9,84 12,43 Altura Líquida (m) 11,87 11,94 11,88 Vazão (m³/s) 0,13962 0,19959 0,24246 Torque no gerador (N.m) 416,39 616,13 740,54 Perda de Torque (N.m) 4,4 4,52 3,85 Velocidade angular do rotor (rpm) 332,84 332,59 332,59 Pressão absoluta na entrada da Voluta (kPa) 218,08 215,57 212,38 Pressão absoluta de saída no tubo de descarga (kPa) 113,17 111,13 109,59 Eficiência Hidráulica (%) 90,13% 92,39% 91,71% Massa específica da água (kg/m³) 999,8 999,8 999,8 Viscosidade Cinemática da água (m²/s) 9,57E-07 9,57E-07 9,57E-07 Aceleração gravitacional (m/s²) 9,82 9,82 9,82 Fonte: NTNU (2021). Para o cálculo da eficiência, é necessário o cálculo da altura líquida disponível (H), sendo determinada pela Equação (8): H = 𝑃1𝑎𝑏𝑠 − 𝑃2𝑎𝑏𝑠 𝜌 𝑔 + 𝑉1 2 − 𝑉2 2 2 𝑔 + 𝑧 (8) onde: P1abs e P2abs são valores absolutos referentes as pressões de entrada na voluta (obtido através dos resultados CFD) e saída no tubo de descarga (dado de entrada para simulação) em [Pa], respectivamente. Os valores de V1 e V2 são as 25 velocidades absolutas na entrada e saída em [m/s], respectivamente, obtidas da relação entre a vazão volumétrica e a seção da área analisada. ρ é a massa específica [kg/m3], g é a aceleração gravitacional em [m/s2]. A diferença de altura entre a entrada da voluta e a saída do tubo de descarga (z) é de 1,075 m (36). O cálculo da eficiência é obtido pela Equação (9) (1): 𝜂 = 𝜔 𝑇 𝜌 �̇�𝑔 𝐻 (9) onde: ω é a rotação do rotor em [rad/s]. T é o torque gerado pelo rotor em [N.m], e é calculado pelo Ansys CFX 19.0 com base no centro de massa dos eixos x, y e z (2). A vazão volumétrica 𝑄 ̇ é utilizada como condição de contorno em [m³/s]. E a altura líquida disponível é obtida pela Equação (8). 3.2 DOMÍNIO COMPUTACIONAL E CONDIÇÕES DE CONTORNO A geometria investigada neste trabalho considera a voluta, com duas palhetas direcionadoras móveis sob a condição de periodicidade (total de 28 pás), uma pá inteira e um splitter (traduzindo do inglês, pá-auxiliar) também sob a condição de periodicidade (total de 15 conjuntos) e o tubo de descarga, conforme mostrado na Figura 18. As posições das geometrias e malhas obtidas pelo Workshop já estavam predefinidas e foram mantidas. A Figura 19 apresenta os domínios computacionais e as interfaces entre os domínios rotativos e estacionários. A Interface I representa a conexão entre a voluta (domínio estacionário) e as palhetas direcionadoras móveis (domínio estacionário). A Interface II representa a conexão entre as palhetas direcionadoras móveis (domínio estacionário) e o rotor (domínio rotativo). Já a Interface III representa a conexão entre o rotor (domínio rotativo) e o tubo de descarga (domínio estacionário). Nas interfaces foram adotadas a condição frozen-rotor (traduzindo em português, rotor parado). Este tipo de interface move o componente conforme a configuração da simulação, porém a conexão entre as duas malhas tem uma posição relativa fixa, produzindo uma solução estacionária para um problema de referência múltipla. 26 Logo, se um domínio de malha deve se mover as transformações das equações são feitas apropriadamente para manter fixa a referência entre elas (2). Figura 18: Domínio computacional. Fonte: Autoria Própria. Figura 19: Domínio computacional e Interfaces. Fonte: Autoria Própria 27 A Tabela 2 apresenta as condições de contorno aplicadas ao domínio computacional da voluta, conforme mostrado na Figura 20. Onde a Interface I apresenta uma distribuição uniforme devido à sua forma construtiva, pois há uma condição de parede na porção mais delgada no final da voluta. Isso limita o escoamento, evitando assim uma possível recirculação. Tabela 2: Condições de contorno no domínio computacional da voluta. Contornos Condições Entrada Vazão mássica prescrita fixa Intensidade turbulenta = 5% Temperatura constante Paredes e palhetas direcionadoras fixas Sem deslizamento Adiabático Parede lisa Interface I Interface GGI (abreviação do inglês General Grid Interface) com palhetas direcionadoras móveis Condição de Frozen Rotor Fonte: Autoria Própria. Figura 20: Domínio computacional da voluta. Fonte: Autoria Própria 28 A Tabela 3 apresenta as condições de contorno no domínio computacional nas palhetas direcionadoras móveis, conforme mostra o conjunto de imagens na Figura 21. Tabela 3: Condições de contorno no domínio computacional das palhetas direcionadoras móveis. Contornos Condições Interface I Interface GGI com voluta Condição de Frozen Rotor Parede e palhetas direcionadoras móveis Sem deslizamento Adiabático Parede lisa Periodicidade Periódico rotacional Fluxo de interface conservativa Interface II Interface GGI com conjunto de pás do rotor Condição de Frozen Rotor Fonte: Autoria Própria. Figura 21: Domínio computacional das palhetas direcionadoras móveis. Fonte: Autoria Própria 29 A Tabela 4 apresenta as condições de contorno no domínio computacional do conjunto de pás do rotor, mostrado pelo conjunto de imagens na Figura 22. Tabela 4: Condições de contorno no domínio computacional do conjunto de pás do rotor. Contornos Condições Superfície da pá, Hub (Cubo) e Shroud (Capa) Sem deslizamento Adiabático Parede lisa Interface II Interface GGI com palhetas direcionadoras móveis Condição de Frozen Rotor Interface III Interface GGI com tubo de descarga Condição de Frozen Rotor Periodicidade Periódico rotacional Fluxo de interface conservativa Fonte: Autoria Própria. Figura 22: Domínio computacional do conjunto de pás do rotor. Fonte: Autoria Própria. 30 A Tabela 5 apresenta as condições de contorno no domínio computacional do tubo de descarga, conforme mostrado pelo conjunto de imagens na Figura 23. Tabela 5: Condições de contorno no domínio computacional do tubo de descarga. Contornos Condições Interface III Interface GGI com conjunto de pás do rotor Condição de Frozen Rotor Cone Com deslizamento Adiabático Saída Pressão Relativa constante Parede Sem deslizamento Adiabático Parede lisa Fonte: Autoria Própria. Figura 23: Domínio computacional do tubo de descarga. Fonte: Autoria Própria. 31 3.3 DISCRETIZAÇÃO DO MODELO COMPUTACIONAL As malhas disponibilizadas pela NTNU, embora compatível com o software Ansys para leitura e utilização em simulações numéricas, não permitem imediata edição devido ao formato do arquivo (.msh). Para a investigação da inserção de dispositivos passivos sobre as superfícies do tubo de descarga e rotor, é necessário adequar as geometrias disponibilizadas pelo Workshop Francis-99 e gerar novas malhas para a realização das simulações numéricas. Assim, é possível inserir dispositivos passivos para verificar seu impacto sobre a eficiência e torque da turbina hidráulica, além do estudo qualitativo do vórtice de corda. 3.3.1 Discretização da geometria do Rotor As malhas do rotor foram geradas considerando dois tipos de elementos, sendo os prismáticos inseridos próximos à parede para melhor caracterizar o desenvolvimento da camada limite dinâmica, e os tetraédricos inseridos no restante do domínio computacional, devido a sua excelente característica de adaptação a geometrias complexas. Para o estudo da densidade de malha, adotou-se um fator de refinamento acima de 1,3, conforme sugere Celik (4). Neste caso, para a modelagem dos elementos prismático, manteve-se constante a altura total próxima a parede das pás, com 8 subdivisões e uma taxa de crescimento entre as camadas de 1,4. Quanto a qualidade da malha, segundo recomendações no Manual do software Ansys, é recomendado uma ortogonalidade maior do que 0,1, sendo que quanto mais próximo da unidade melhor será a qualidade da malha (2). O parâmetro de ortogonalidade entre células é obtido a partir do vetor normal da face (Vetores A), Figura 24, o vetor do centroide da célula ao centroide das células adjacentes (Vetores C), e o vetor do centroide da célula à cada centroide das faces da célula (Vetores f). 32 Figura 24: Esquema apresentando a definição de ortogonalidade. Fonte: Ansys (2011). A qualidade ortogonal para cada célula é então calculada através dos vetores 𝐴𝑖 ⃗⃗ ⃗, 𝑓𝑖⃗⃗ e 𝐶𝑖 ⃗⃗ ⃗, considerando o valor máximo entre os parâmetros Af e AC (max (Af, AC)) definidos nas Equações (10) e (11), respectivamente, a seguir: 𝐴𝑓 = 𝐴𝑖 ⃗⃗ ⃗ . 𝑓𝑖⃗⃗ |𝐴𝑖 ⃗⃗ ⃗| |𝑓𝑖⃗⃗ | (10) 𝐴𝐶 = 𝐴𝑖 ⃗⃗ ⃗ . 𝐶𝑖 ⃗⃗ ⃗ |𝐴𝑖 ⃗⃗ ⃗| |𝐶𝑖 ⃗⃗ ⃗| (11) Com os valores de ortogonalidade de cada célula obtidos pelos parâmetros Af e AC, o software Ansys Meshing utiliza o valor máximo de cada célula e realiza uma média aritmética simples, sendo possível ter um parâmetro numérico de avaliação (2). A Tabela 6 apresenta as informações sobre os refinos da malha do tubo de descarga, e a respectiva ortogonalidade mínima e média das malhas. Por ser uma geometria complexa com superfícies irregulares, os índices de ortogonalidade da malha na região do rotor podem ser considerados aceitáveis para realizar as simulações. 33 Tabela 6: Características da discretização do rotor. Malha Número de Elementos Ortogonalidade mínima Ortogonalidade média N1 (Refinada) 7.525.640 0,11 0,78 N2 (Intermediária) 3.306.725 0,10 0,76 N3 (Grosseira) 1.493.426 0,10 0,73 Fonte: Autoria Própria. A Figura 25 apresenta as malhas geradas para o canal do rotor (a-c), e uma vista em corte (d-f) com detalhes das camadas de prismas sobre a superfície da pá. Figura 25: Refino de malha do canal do rotor. Fonte: Autoria Própria. 34 3.3.2 Discretização da geometria do Tubo de Descarga As malhas no tubo de descarga foram geradas com três tipos de elementos, sendo elementos tetraédricos na região de entrada do cone devido à adequada adaptação à geometrias complexas, elementos prismáticos próximo a parede para caracterização do desenvolvimento da camada limite dinâmica, e no restante do tubo definiu-se elementos hexaédricos, reduzindo a quantidade de elementos da malha sem prejuízo à solução da modelagem e ao processo de convergência. Para a realização do refino de malha, novamente foi definido um fator de refinamento maior do que 1,3, seguindo recomendações de Celik (4), sendo definido para os elementos prismáticos 8 camadas próximo a parede, com uma taxa de crescimento entre camadas de 1,4 e mantida a altura do primeiro elemento em 0,15 mm. A Tabela 7 apresenta as informações sobre o refino de malha no tubo de descarga e as respectivas ortogonalidades mínimas e médias das malhas. Como pode ser observado, a ortogonalidade média está próxima de 1 e a mínima maior do que 0,1, indicando que as malhas estão adequadas segundo as recomendações em (2). Tabela 7: Propriedades do refino de malha do tubo de descarga. Malha Número de Elementos Ortogonalidade mínima Ortogonalidade média N1 (Refinada) 8.172.750 0,19 0,93 N2 (Intermediária) 3.548.905 0,29 0,93 N3 (Grosseira) 1.569.722 0,23 0,94 Fonte: Autoria Própria. A Figura 26 mostra as malhas geradas em corte transversal na parte do cone (a-c) e uma vista superior da malha hexaédrica (d-f). 35 Figura 26: Refino de malha do tubo de descarga. Fonte: Autoria Própria. 3.4 ANÁLISE PRELIMINAR O objetivo desta análise preliminar é verificar se a análise de convergência de malha, através da abordagem proposta por Celik (4), poderia ser realizada inicialmente por meio de uma modelagem estacionária em contrapartida a uma abordagem transitória. O torque e a eficiência da turbina hidráulica Francis foram avaliados e comparados. Além disso, uma análise qualitativa do vórtice de corda também é realizada. Para este estudo, foi utilizada a malha padrão disponibilizada pelo Workshop Francis-99, sendo que para a abordagem transitória a análise é realizada para a décima revolução. Para a configuração do regime transitório foram utilizadas as configurações de simulação numéricas adotadas por Tran (44), onde a turbina Francis-99 foi simulada sob a condição operacional PL, com o objetivo de identificar o vórtice de corda, onde o tempo total avaliado é de 1,8 segundo, equivalente a dez revoluções da turbina em condições reais, com passo de tempo de 5E-4 segundo 36 (equivalente ao avanço angular de 1°), com critério de convergência de erro residual RMS (sigla derivada do inglês Root-Mean-Square, traduzindo Raiz Quadrática Média) de 1E-5. A solução transitória foi inicializada a partir da solução em regime estacionário. Conforme pode ser verificado na Tabela 8, uma pequena diferença percentual é observada entre os valores de torque e eficiência para os regimes estacionário e transitório para os três pontos de operação. Tabela 8: Comparação entre hipóteses de regimes estacionário e transitório. Condição Regime Torque (N.m) Diferença Torque (%) Eficiência (%) Diferença Eficiência (%) PL Estacionário 475,90 0,83 87,11 0,09 Transitório 479,84 87,03 BEP Estacionário 677,59 1,38 87,89 0,04 Transitório 686,97 87,93 HL Estacionário 803,43 1,48 86,60 0,18 Transitório 815,30 86,76 Fonte: Autoria Própria. A partir dos dados apresentados na Tabela 8, sob a condição operacional PL, uma análise qualitativa do vórtice de corda foi realizada utilizando o método Q- criterion com valor fixo de Q = 300 s-2, inserindo este valor em uma função de pós- processamento no software CFD-Post do Ansys 19.0, atribuindo a grandeza de velocidade. Como pode ser visto na Figura 27, o vórtice de corda não foi completamente caracterizado com a malha disponibilizada pelo Workshop Francis- 99, conforme reportado por Tran (44). Desta forma, para a caracterização do vórtice de corda tanto para o regime estacionário quanto para o regime transitório, foi desenvolvida uma malha do rotor mais refinada quando comparada à malha padrão disponibilizada pelo Workshop Francis-99. A malha analisada nesta avaliação preliminar é semelhante à malha intermediária detalhada na seção 3.3.1. Como pode ser observado na Figura 28, com uma malha mais refinada o vórtice de corda é melhor caracterizado, semelhante ao mostrado na Figura 10. Além disso, embora o vórtice de corda apresente diferença qualitativa em seu 37 formato considerando uma abordagem estacionária comparada com a abordagem transitória, sendo ambas com domínios periódicos, é razoável a estratégia indicada na seção 3, onde propõem-se iniciar o estudo da inserção de dispositivos passivos sobre as pás do rotor e tubo de descarga a partir de uma abordagem estacionária e posteriormente submeter a solução promissora a uma análise transitória. Figura 27: Malha padrão na condição de operação PL com Q-criterion 300 s-2. Fonte: Autoria Própria. 38 Figura 28: Comparação Estacionário e Transitório no Regime PL. Fonte: Autoria Própria. O número de Courant, é um parâmetro que foi avaliado para a abordagem transitória, que verifica a relação do passo de tempo, refino de malha e a velocidade local do escoamento, através da equação (2): 𝐶𝑜𝑢𝑟𝑎𝑛𝑡 = 𝑢 Δ𝑡 Δ𝑥 (12) onde: u é a velocidade local do escoamento em [m/s], Δ𝑡 é o passo de tempo da simulação transitória em [s], e Δ𝑥 é o tamanho do elemento da malha local em [m]. Na documentação do software Ansys, é recomendado que o passo de tempo seja definido para que resulte em um número de Courant suficientemente pequeno. Entretanto, não é recomendado um valor ideal para o número de Courant e nem o local a ser avaliado, pois, naturalmente, depende de cada modelagem computacional (2). Contudo, considerando outras temáticas considerando domínios rotativos e estacionários, decidiu-se avaliar o número de Courant nas interfaces dos domínios 39 computacionais, conforme a Figura 29. A Figura 29 (a) apresenta o número de Courant na Interface I da voluta. A Figura 29 (b) apresenta a análise da Interface I das palhetas direcionadoras móveis, enquanto a Figura 29 (c) apresenta a Interface II do mesmo domínio. A Figura 29 (d) apresenta o número de Courant da Interface II do rotor, enquanto a Figura 29 (e) apresenta a análise da Interface III do rotor. Já a Figura 29 (f) apresenta a análise do número de Courant da Interface III do tubo de descarga. Figura 29: Número de Courant nas interfaces dos domínios computacionais. Fonte: Autoria Própria. 40 3.5 ESTUDO DA DENSIDADE DE MALHA O estudo quanto à densidade de malha foi realizado através da abordagem proposta por Celik (4), conhecido como Grid Convergence Index (GCI). Este método é baseado nas Extrapolações de Richardson, sendo que para seu desenvolvimento é necessário o uso de três diferentes refinos de malhas, conforme detalhado na seção 3.3.1 e 3.3.2. A Tabela 9 apresenta as informações sobre o número de elementos total considerando a somatória dos elementos no rotor e no tubo de descarga, além do respectivo valor do fator de refinamento, conforme recomendação de Celik (4). Tabela 9: Informações das malhas para o GCI. Malha Elementos Fator de Refinamento N1 15.698.390 1,32 N2 6.855.630 1,31 N3 3.063.148 * Fonte: Autoria Própria. O critério de convergência adotado para todas as variáveis, pressão, velocidade e modelo de turbulência SST k- ω, é o erro residual médio (RMS) de 1E- 5. Porém, o processo de convergência para as condições operacionais PL e HL apresentaram-se ser mais instáveis quando comparados ao BEP, característica também reportada por Stoessel (42). Desta forma, com o objetivo de reduzir o tempo de processamento, avaliou-se adicionalmente, ao longo do processo de convergência da solução, a estabilidade dos principais parâmetros de saída (Torque e Eficiência). Conforme mostrado na Figura 30 e Figura 31, para a condição operacional PL, e considerando a malha intermediária (N2), apenas a pressão (indicado por RMS P-Mass) superou o erro de 1E-5. E conforme mostrado nas Figura 30, Figura 31, Figura 32 e Figura 33 a partir da iteração 1500 os resultados oscilam em uma estreita faixa, indicando razoável estabilidade, permitindo o truncamento da solução numérica. 41 Figura 30: Critério de convergência do erro residual médio (RMS) da pressão e velocidade. Fonte: Autoria Própria. Figura 31: Critério de convergência do erro residual médio (RMS) de SST k-ω. Fonte: Autoria Própria. 1,0E-06 1,0E-05 1,0E-04 1,0E-03 1,0E-02 1,0E-01 1,0E+00 0 500 1000 1500 2000 E rr o R e s id u a l (R M S ) Iteração RMS P-Mass RMS U-Mom RMS V-Mom RMS W-Mom 1,0E-06 1,0E-05 1,0E-04 1,0E-03 1,0E-02 1,0E-01 1,0E+00 0 500 1000 1500 2000 E rr o R e s id u a l (R M S ) Iteração RMS Turbulência (Epsilon) RMS Turbulência (Ômega) 42 Figura 32: Perfil de convergência da eficiência. Fonte: Autoria Própria. Figura 33: Perfil de convergência do torque. Fonte: Autoria Própria. 84,0 84,5 85,0 85,5 86,0 86,5 87,0 87,5 88,0 0 500 1000 1500 2000 E fi c iê n c ia ( % ) Iteração Eficiência 455,0 455,5 456,0 456,5 457,0 457,5 458,0 458,5 459,0 459,5 460,0 0 500 1000 1500 2000 T o rq u e ( N .m ) Iteração Torque 43 O formato do vórtice de corda também foi avaliado qualitativamente a partir do método Q-criterion com valor de Q = 300 s-2, sob a condição operacional PL e refino de malha intermediária (N2), conforme mostrado na Figura 34. Como pode ser observado, não houve diferença significativa no formato do vórtice de corda para as iterações 1.500 e 2.000, permitindo concluir que há estabilidade no seu formato a partir da iteração 1500. Para as demais condições operacionais (BEP e HL), conclusões semelhantes também foram observadas. Figura 34: Critério de Convergência por análise qualitativa. Fonte: Autoria Própria. Desta forma, o cálculo da métrica GCI foi realizado considerando as informações de malha mostradas na Tabela 9 e Tabela 10. Como pode ser visto, os valores de GCI21 (incerteza numérica da malha mais refinada) e GCI32 (incerteza 44 numérica da malha intermediária) são muito semelhantes, indicando que o uso da malha intermediária não traz prejuízo aos resultados de torque e eficiência da turbina hidráulica. O cálculo de GCI foi realizado considerando os valores apresentados na Tabela 11. Tabela 10: Análise de incerteza devido a densidade de malha para o torque e eficiência. Operação GCI Torque (%) Eficiência (%) PL GCI21 0,003226% 0,121825% GCI32 0,000280% 0,150832% BEP GCI21 0,000005% 0,002397% GCI32 0,000584% 0,000418% HL GCI21 0,000932% 0,026868% GCI32 0,000059% 0,079906% Fonte: Autoria Própria. Tabela 11: Comparação entre os resultados a partir das malhas N1, N2 e N3 e Workshop Francis-99. PL BEP HL Malha Torque (N.m) Eficiência (%) Torque (N.m) Eficiência (%) Torque (N.m) Eficiência (%) N1 456,944 86,290 648,300 87,082 768,660 86,255 N2 456,820 86,270 648,303 87,074 768,745 86,220 N3 456,810 86,246 647,973 87,073 768,740 86,320 Fonte: Autoria Própria. 3.6 VALIDAÇÃO NUMÉRICA A validação da abordagem numérica é feita pelo confronto com os dados experimentais a partir do Workshop Francis-99. Além disso, os resultados numéricos do presente trabalho foram comparados os resultados das simulações numéricas desenvolvidas por Minakov (32) e Gravilov (12), os quais utilizaram o modelo de turbulência k-ω SST. Além disso, Minakov (32) modelou a turbina completa, enquanto Gravilov (12) não utilizou a voluta, inserindo perfis de velocidades na entrada das palhetas direcionadoras móveis. 45 A Figura 35 mostra que os resultados do torque para o presente trabalho indicam similaridade quando comparada aos dados experimentais e também em relação aos trabalhos de Minakov (32) e Gravilov (12), embora diferentes abordagens de modelagem tenham sido usadas. A análise das eficiências obtidas neste trabalho também são comparadas aos dados experimentais a partir do Workshop Francis-99 e também ao trabalho de Minakov (32). Como pode ser visto na Figura 36, os resultados obtidos no presente trabalho indicam valores próximos aos experimentais e comportamento parecido aos obtidos por Minakov (32), embora existam diferenças nas condições de contorno e domínio computacional que possam justificar os perfis. Entretanto, para os propósitos deste trabalho, a modelagem pode ser considerada satisfatória. Figura 35: Comparação de resultados numéricos e experimental para o Torque (T). Fonte: Autoria Própria 46 Figura 36: Comparação de resultados numéricos e experimental para a Eficiência. Fonte: Autoria Própria 47 4 RESULTADOS E DISCUSSÕES Inicialmente, alguns tipos de dispositivos passivos inseridos sobre as superfícies do tubo de descarga e também sobre as pás foram avaliados quanto ao impacto na formação do vórtice de corda, considerando a condição operacional PL em regime estacionário. 4.1 RESULTADOS PARA DISPOSITIVOS PASSIVOS INSTALADOS NO TUBO DE DESCARGA Para avaliar a possibilidade de mitigação do vórtice de corda foram inseridos dispositivos passivos arbitrariamente no tubo de descarga com o objetivo de direcionar o escoamento do fluido para a região central do tubo de descarga, perturbando a zona de estagnação, conforme conjunto 1 de possíveis dispositivos passivos mostrados na Figura 37. Nas alternativas mostradas pelo conjunto 2, Figura 38, o objetivo é semelhante ao do conjunto 1, com geometrias em locais diferentes, embora o dispositivo passivo DT_5 se encontre fixo sobre o cone. 48 Figura 37: Dispositivos Passivos no tubo de descarga (conjunto 1). Fonte: Autoria Própria. 49 Figura 38: Dispositivos Passivos no tubo de descarga (conjunto 2). Fonte: Autoria Própria. As nomenclaturas das geometrias são indicadas abaixo: • Original: Sem inserção de dispositivo passivo. • DT: abreviação do inglês Draft Tube, traduzindo Tubo de Descarga. Conforme pode ser observado na Tabela 12, todas as alternativas apresentaram alteração nos valores de eficiência e torque da turbina hidráulica na condição PL, com destaque para o dispositivo passivo DT_3, o qual reduziu o valor da eficiência de forma mais significativa. Os demais dispositivos passivos mantiveram praticamente os mesmos valores de torque e eficiência, com diferenças praticamente desprezíveis. 50 Tabela 12: Resultados das simulações com DP no tubo de descarga. Geometria Torque (N.m) Eficiência (%) Original 456,82 * 86,270 * DT_1 456,80 -0,004% 86,280 0,012% DT_ 2 456,81 -0,002% 86,243 -0,031% DT_ 3 456,80 -0,004% 85,852 -0,485% DT_ 4 456,80 -0,004% 86,216 -0,063% DT_ 5 456,83 0,002% 86,274 0,005% DT_ 6 456,82 0,000% 86,274 0,005% Fonte: Autoria Própria. A Figura 39 e Figura 40 mostram uma análise qualitativa do vórtice de corda para cada alternativa avaliada nesta seção. Analisando os resultados para o conjunto 1, mostrados na Figura 39, percebe-se que o dispositivo passivo DT_3 gera um vórtice de corda muito intenso após o estreitamento causado pela geometria passiva, reduzindo em aproximadamente 0,5% a eficiência da turbina, indicando que esta alternativa não é promissora para as análises futuras. Os demais dispositivos passivos não impactaram significativamente a formação do vórtice de corda, provocando alterações na parede do tubo de descarga. 51 Figura 39: Q-criterion de 300 s-2 dos dispositivos passivos instalados no tubo de descarga (conjunto 1). Fonte: Autoria Própria. A análise do conjunto 2, Figura 40, indica que o dispositivo passivo DT_4 promove significativa alteração no escoamento próximo a parede do tubo de descarga, fragmentando o vórtice de corda na região a jusante, tornando esse dispositivo passivo promissor para um melhor detalhamento futuro. Os demais dispositivos passivos do conjunto 2 não impactaram significativamente o vórtice de corda, requerendo também mais análises com diferentes tamanhos e posições destas geometrias para a tomada de decisão mais contundente. 52 Figura 40: Q-criterion de 300 s-2 dos dispositivos passivos instalados no tubo de descarga (conjunto 2). Fonte: Autoria Própria. Contudo, a inserção de dispositivos passivos no tubo de descarga não gerou impacto significativo no torque do rotor e na eficiência, com exceção do dispositivo passivo DT_3. O vórtice de corda se mostra mais suscetível às mudanças com a inserção de geometrias maiores e próximos da zona de estagnação. Entretanto, nenhuma alternativa com dispositivos passivos no tubo de descarga se mostrou promissora. 53 4.2 RESULTADOS PARA DISPOSITIVOS PASSIVOS SOBRE O ROTOR TIPO 1 Uma alternativa inicial para avaliar a possibilidade de mitigar o vórtice de corda com a inserção de dispositivos passivos arbitrariamente sobre as superfícies das pás do rotor é apresentada na Figura 41, com o objetivo de investigar a dinâmica do escoamento sobre a forma do vórtice de corda conduzindo o escoamento do fluido para a região central na saída do rotor, com a expectativa de reduzir a formação da zona de estagnação no tubo de descarga. Figura 41: Dispositivos Passivos (DP) Tipo 1, do tipo barreira no rotor. Fonte: Autoria Própria. Os dispositivos foram posicionados próximos ao Hub e Shroud, conforme nomenclatura a seguir: • Original: Sem inserção de barreiras. • RU: abreviação do inglês Runner, traduzindo Rotor. • SH: abreviação do inglês Shroud, traduzindo Capa. • HUB: traduzindo do inglês Cubo. O dispositivo passivo RU_SH_2 é o único com perfil retangular com altura menor em comparação aos demais dispositivos passivos que possuem perfil em meia circunferência, conforme mostrado na Figura 42. 54 Figura 42: Perfis dos dispositivos passivos do Tipo 1 sobre o rotor. Fonte: Autoria Própria. Conforme pode ser observado na Tabela 13, todas as alternativas de dispositivos passivos do tipo 1 instalados sobre o rotor resultaram em modificações nos valores de eficiência e/ou torque da turbina hidráulica na condição PL. O dispositivo passivo RU_SH_1 apresentou maior redução do torque, e a configuração RU_HUB_SH_1 apresentou a maior redução na eficiência. Assim, as alternativas preliminarmente estudadas não demonstraram adequada perspectiva quanto ao impacto no torque e eficiência da turbina hidráulica. Tabela 13: Resultados das simulações com DP do tipo barreira no rotor. Geometria Torque (N.m) Eficiência (%) Original 456,82 * 86,270 * RU_SH_1 453,02 -0,832% 85,690 -0,673% RU_SH_2 455,65 -0,256% 85,943 -0,379% RU_HUB_1 456,82 0,000% 86,114 -0,181% RU_HUB_SH_1 453,67 -0,690% 85,565 -0,817% Fonte: Autoria Própria. Além disso, uma análise qualitativa do vórtice de corda foi realizada considerando cada alternativa de dispositivo passivo do tipo 1 instalada sobre o rotor, considerando o método do Q-criterion com valor de Q = 300 s-2, conforme mostrado na Figura 43. Nota-se que o dispositivo passivo RU_SH_1 é aquele que mais impacta no formato do vórtice de corda, aumentando seu tamanho com relação ao modelo original. O dispositivo passivo RU_HUB_1 foi a alternativa que 55 apresentou uma pequena redução no vórtice de corda e na dinâmica do escoamento próximo à parede do tubo de descarga. Figura 43: Q-criterion de 300 s-2 dos DP’s do tipo barreira no rotor. Fonte: Autoria Própria. Verificou-se que a inserção deste tipo 1 de dispositivos passivos instalados sobre o rotor impactou na redução da eficiência e no torque da turbina hidráulica. Além disso, o vórtice de corda se torna suscetível à mudança de comportamento com a inserção dos dispositivos passivos, embora não necessariamente para a sua mitigação. 56 4.3 RESULTADOS PARA OS DISPOSITIVOS PASSIVOS SOBRE O ROTOR TIPO 2 Outra abordagem de dispositivo passivo instalado sobre o rotor é mostrado na Figura 44, onde o objetivo continua sendo a modificação da dinâmica do escoamento, centralizando a saída da água do rotor na região de estagnação, minimizando a formação do vórtice de corda. A Figura 44 mostra o rotor secionado transversalmente para facilitar a visualização do dispositivo, destacado em vermelho. O formato, tamanho e posição desta alternativa foram inicialmente modelados para permitir avaliar o comportamento do vórtice de corda. Sua construção foi desenvolvida para que não seja permitido o escoamento no seu interior, mas apenas na região externa no canal da turbina hidráulica. Figura 44: Detalhe Dispositivo Passivo Tipo 2, do tipo barreira no rotor. Fonte: Autoria Própria. Como pode ser visto na Figura 45, o vórtice de corda na condição de operação PL permaneceu centralizado, com resultados satisfatórios quando comparados às possibilidades abordadas nas seções 4.1 e 4.2. Porém a velocidade de escoamento aumentou significativamente devido ao estreitamente na área de passagem no rotor. 57 Figura 45: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2. Fonte: Autoria Própria. Considerando a análise preliminar apresentada no item 3.4 deste trabalho, o dispositivo passivo do Tipo 2 instalado sobre o rotor se torna interessante no ponto de vista de mitigação do vórtice de corda, embora características geométricas ainda necessitem ser avaliadas, além do seu impacto sobre as condições de contorno BEP e HL. 58 4.4 RESULTADOS PARA O DISPOSITIVO PASSIVO TIPO 2 COM VARIAÇÃO NOS PARÂMETROS GEOMÉTRICOS Nesta análise, os dispositivos passivos foram instalados sobre a superfície da pá do rotor próximo ao bordo de fuga, considerando diferentes tamanhos e posição e para as três condições de operação: PL, BEP e HL. Conforme discutido na seção 3.4, inicialmente um regime estacionário é avaliado, e a abordagem transitória será discutida apenas para a configuração mais promissora para a mitigação do vórtice de corda, que também apresente pequeno impacto sobre o torque e a eficiência, embora seja desejável o aumento destes parâmetros. A Figura 46 mostra o rotor em corte com detalhe em vermelho do dispositivo passivo do tipo 2. A água, ao adentrar o rotor, escoa por entre as pás principais e auxiliares, onde o dispositivo passivo desvia o escoamento para a região central do rotor. As variáveis investigadas, como mostrado na Figura 46, são: comprimento (L) e ângulo (θ). A altura (h) foi mantida constante, conforme valores mostrados na Tabela 14. Figura 46: Detalhe Dispositivo Passivo Tipo 2, do tipo barreira no rotor com parâmetros. Fonte: Autoria Própria. 59 Tabela 14: Parâmetros de variação dos dispositivos passivos. Comprimento (L) Ângulo (θ) RU_25_60 25 mm 60° RU_25_45 25 mm 45° RU_25_30 25 mm 30° RU_25_15 25 mm 15° RU_12,5_60 12,5 mm 60° RU_12,5_45 12,5 mm 45° RU_12,5_30 12,5 mm 30° RU_12,5_15 12,5 mm 15° Fonte: Autoria Própria. O ângulo 𝜃 representa a inclinação do dispositivo passivo em relação a um plano paralelo à saída do rotor, ou seja, quanto menor o ângulo, maior será a restrição no canal. A dinâmica do escoamento devido à configuração de dispositivo passivo tipo 2 com parâmetros sobre o rotor foi inicialmente feita a partir da comparação com uma configuração sem dispositivo passivo e outra com um dispositivo passivo (RU_25_15) que apresenta a maior redução de área da passagem de água, conforme mostrado na Figura 47. É possível observar que o perfil de velocidade, plotado em um plano transversal na saída do rotor, é significativamente alterado devido ao dispositivo passivo, com um aumento de velocidade próximo ao centro do rotor (cubo). Figura 47: Comparação do perfil de velocidade com e sem barreira. Fonte: Autoria Própria. 60 A análise segue considerando as configurações de dispositivos passivos mostrados na Tabela 14, separados por condição operacional e comprimento do dispositivo passivo. A configuração sem dispositivo passivo é nomeada como “Original”. Para a condição operacional PL e dispositivo passivo de comprimento de 25 mm, a Figura 48, mostra os vórtices de corda formados entre o rotor e o tubo de descarga em vista lateral e também em vista superior, com a mesma escala de cor da velocidade do escoamento com um Q-criterion de 300 s-2. Pode-se observar que todas as configurações apresentaram a formação do vórtice de corda centralizado e redução da dinâmica do escoamento na parede, onde o vórtice de corda mais curto é verificado para o dispositivo passivo RU_25_60, enquanto a configuração RU_25_15 gera um vórtice de corda mais delgado. Figura 48: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=25 mm, na condição de operação PL. Fonte: Autoria Própria. Utilizando a mesma abordagem de análise, na Figura 49 são apresentados os resultados para os dispositivos passivos de comprimento 12,5 mm e regime de 61 operação PL. Observa-se que a única configuração que apresentou um vórtice de corda centralizado e curto foi a RU_12,5_15, embora a dinâmica do escoamento na parede aumentou em todas as configurações. A configuração RU_12,5_60 e o RU_12,5_45 aumentaram o comprimento do vórtice de corda, enquanto o RU_12,5_30 reduziu sutilmente o comprimento do vórtice de corda, ficando mais centralizado e delgado. Figura 49: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=12,5 mm, na condição de operação PL. Fonte: Autoria Própria. Para o regime de operação BEP, utilizando os dispositivos passivos de comprimento de 25 mm, é possível observar na Figura 50 que o dispositivo passivo RU_25_60 apresenta uma mudança sutil apenas sobre o vórtice na parede, similar à estrutura gerada pela configuração Original, enquanto a configuração RU_25_45 apresentara mudança sutil sobre o vórtice de corda. Para a configuração RU_25_30, um vórtice de corda mais curto é verificado, porém centralizado e robusto. O dispositivo passivo RU_12_15 apresentou o vórtice de corda mais curto e robusto, com o surgimento de estruturas de vórtices próximo à parede. 62 Figura 50: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=25 mm, na condição de operação BEP. Fonte: Autoria Própria. Na Figura 51, regime de operação BEP, pode ser observado que todos os dispositivos passivos de comprimento 12,5 mm mantiveram o padrão do vórtice de corda central que são similares ao gerado pela configuração Original, embora tenha se observado o aumento da dinâmica do escoamento próximo à parede quando comparados à configuração Original. 63 Figura 51: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=12,5 mm, na condição de operação BEP. Fonte: Autoria Própria. Para a condição de operação HL e dispositivos passivos de comprimento 25 mm, pode-se observar na Figura 52 que o dispositivo passivo RU_25_60 gerou um vórtice de corda central similar à configuração Original, com alteração sutil na dinâmica do escoamento próximo à parede. As demais variações de ângulo de ângulo deste dispositivo passivo de 25 mm geraram vórtices de corda central mais robusto e também com vórtices próximos à parede. 64 Figura 52: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=25 mm, na condição de operação HL. Fonte: Autoria Própria. Para o mesmo regime de operação HL, porém com dispositivo passivo de comprimento de 12,5 mm, a Figura 53 mostra que os dispositivos passivos RU_12,5_60 e RU_12,5_45 apresentaram mudanças sutis em relação à configuração Original. As demais configurações apresentaram um vórtice de corda central mais robusto. 65 Figura 53: Q-criterion de 300 s-2 do Dispositivo Passivo Tipo 2 com L=12,5 mm, na condição de operação HL. Fonte: Autoria Própria. A Tabela 15 mostra os resultados para o torque e eficiência para cada dispositivo passivo e condição de operação. A diferença percentual reportada é em relação aos dispositivos passivos e à geometria Original. Além disso, é destacado na cor verde os melhores torques, e na cor laranja as melhores eficiências. Como observado, os dispositivos passivos que apresentaram os melhores resultados para torque são aqueles que possuem um ângulo de inclinação mais agudo (15º), embora tenham também apresentado os piores resultados para a eficiência. Por outro lado, os dispositivos passivos que apresentaram os melhores resultados para a eficiência são aqueles com ângulos mais abertos (60°). 66 Tabela 15: Resultados de Torque e Eficiência dos dispositivos passivos do tipo 2. Fonte: Autoria Própria. A partir dos resultados anteriores, pode-se inferir qual o dispositivo passivo seria o mais adequado para mitigar o vórtice de corda. Entretanto, deve-se ponderar que o dispositivo passivo não deve impactar negativamente na eficiência e torque, independente da condição operacional, desta forma foram feitas análises separadas para PL, BEP e HL, e consequentemente escolhido o melhor dispositivo passivo dentre os propostos. Partindo de uma análise qualitativa do vórtice de corda na condição de operação PL, conforme mostrado na Figura 48 e Figura 49, os dispositivos passivos RU_25_60 e RU_12,5_15 se destacam por apresentarem o vórtice de corda menor e mais centralizado. Porém, conforme mostrado na Tabela 15 na condição PL, o RU_25_60 apresenta um acréscimo na eficiência, enquanto o RU_12,5_60 há uma perda. Para a condição de operação BEP, conforme apresentado na Figura 50 e na Figura 51, os dispositivos passivos que menos produzem impacto com relação à configuração Original são o RU_25_60 e o RU_12,5_60, embora a Tabela 15 mostre que o dispositivo passivo RU_25_45 se destaca por possuir o maior aumento na eficiência entre aqueles com o mesmo comprimento. Já o dispositivo passivo RU_12,5_60 se destaca por ter uma menor redução na eficiência, também entre os dispositivos de mesmo comprimento. Na condição de operação HL, conforme mostrado na Figura 52 e Figura 53, pode-se perceber que os dispositivos passivos RU_25_60 e RU_12,5_60 obtiveram resultados parecidos em relação à configuração Original. Os resultados mostrados na Tabela 15 também evidenciam esta similaridade, onde são observadas as menores diferenças entre as eficiências. GEOMETRIA Torque (N.m.) Dif. Torque (%) Eficiência (%) Dif. Eficiência (%) Torque (N.m.) Dif. Torque (%) Eficiência (%) Dif. Eficiência (%) Torque (N.m.) Dif. Torque (%) Eficiência (%) Dif. Eficiência (%) Original 456,820 - 86,270 - 648,303 - 87,074 - 768,745 - 86,220 - RU_25_60 460,720 0,854 86,537 0,309 653,070 0,735 86,880 -0,223 779,640 1,417 85,794 -0,494 RU_25_45 466,510 2,121 86,681 0,476 662,040 2,119 86,349 -0,833 791,080 2,905 84,783 -1,666 RU_25_30 470,550 3,006 86,535 0,307 668,310 3,086 85,543 -1,758 799,250 3,968 83,630 -3,004 RU_25_15 472,820 3,502 86,211 -0,068 670,600 3,439 84,680 -2,749 801,830 4,304 82,207 -4,654 RU_12,5_60 455,970 -0,186 86,178 -0,107 645,750 -0,394 87,038 -0,041 769,830 0,141 86,233 0,015 RU_12,5_45 454,860 -0,429 85,948 -0,373 643,910 -0,678 86,808 -0,306 767,350 -0,181 85,938 -0,327 RU_12,5_30 455,460 -0,298 85,855 -0,481 644,240 -0,627 86,427 -0,743 767,580 -0,152 85,325 -1,038 RU_12,5_15 457,540 0,158 85,851 -0,486 647,010 -0,199 86,095 -1,124 770,340 0,207 84,737 -1,720 PL BEP HL 67 Entre todos os dispositivos passivos do tipo 2, a configuração RU_25_60 foi a mais promissora, com destaque para o aspecto qualitativo do vórtice de corda e na eficiência, pois esta geometria reduziu e centralizou o vórtice de corda, além de apresentar menor intensidade turbulenta próximo às paredes para a condição de operação PL quando comparados com o RU_12,5_15. Além disso, o dispositivo passivo RU_25_60 aumentou sutilmente a eficiência da turbina sob a mesma condição de operação. Nas demais condições BEP e HL, este dispositivo gera um comportamento do vórtice de corda similar ao original, porém com uma pequena redução na eficiência, ainda menor se comparado aos demais dispositivos passivos de mesmo comprimento. Com relação ao torque, é observado um aumento nas três condições de operação em relação à original. Desta forma, será utilizado o dispositivo passivo RU_25_60 para avaliação do comportamento do vórtice de corda em regime transitório. 4.5 RESULTADOS DO RU_25_60 NO REGIME TRANSITÓRIO As configurações do regime transitório são as mesmas utilizadas na Seção 3.4 deste trabalho. A análise qualitativa do vórtice de corda considerando uma abordagem em regime transitório para a condição de operação PL é apresentada na Figura 54. Foi definido o ângulo de 90° de rotação do rotor para a visualização dos resulta