Leandro de Souza Alencar

Simulação Numérica de Reinício do

Escoamento de Fluidos Tixotrópicos

DISSERTAÇÃO DE MESTRADO

Dissertação apresentada ao Programa de Pós-

Graduação em Engenharia Mecânica do

Departamento de Engenharia Mecânica do Centro

Técnico Científico da PUC-Rio, como requisito parcial

para o diploma de Mestrado.

Orientadora: Profa. Mônica Naccache

Rio de Janeiro

Fevereiro de 2016

Leandro de Souza Alencar

Simulação Numérica de Reinício do

Escoamento de Fluidos Tixotrópicos

Dissertação apresentada ao Programa de Pós-Graduação em

Engenharia Mecânica do Departamento de Engenharia

Mecânica do Centro Técnico Científico da PUC-Rio, como

requisito parcial para o diploma de Mestrado.

Profa. Mônica Naccache Orientadora

Departamento de Engenharia Mecânica – PUC-Rio

Prof. Paulo Roberto de Souza Mendes Departamento de Engenharia Mecânica – PUC-Rio

Dr. Rafael Mendes Centro de Pesquisa e Desenvolvimento – PETROBRAS

Prof. Marcio da Silveira Carvalho

Coordenador do Centro Técnico Científico da PUC-Rio

Rio de Janeiro

Fevereiro de 2016

Todos os direitos reservados. É proibida a

reprodução total ou parcial do trabalho sem

autorização da universidade, do autor e do

orientador.

Leandro de Souza Alencar Graduou-se em Engenharia de Controle e

Automação pela Universidade Federal de Minas

Gerais (UFMG) em 2008.

Ficha Catalográfica

CDD: 621

Alencar, Leandro de Souza Simulação numérica de reinício do escoamento de fluidos tixotrópicos / Leandro de Souza Alencar; orientadora: Mônica Naccache. – 2016. 80 f. ; 30 cm Dissertação (mestrado)–Pontifícia Universidade Católica do Rio de Janeiro, Departamento de Engenharia Mecânica, 2016. Inclui bibliografia 1. Engenharia mecânica – Teses. 2. Tixotropia. 3. Simulação numérica. 4. CFD. 5. Repartida de escoamento. I. Naccache, Mônica. II. Pontifícia Universidade Católica do Rio de Janeiro. Departamento de Engenharia Mecânica. III. Título.

Agradecimentos

A minha orientadora Mônica Nacchace pelas idéias, discussões e disponibilidade

essenciais para o desenvolvimento deste trabalho.

Ao meu ex-gerente Guilherme Peixoto por incentivar e apoiar minha inscrição ao

programa de mestrado e ao meu gerente Cezar Paulo e coordenador José Ricardo

Montesanti pela compreensão e por permitir minha liberação parcial aos estudos.

A PETROBRAS pelo apoio financeiro.

Ao meu colega Gustavo Moisés por compartilhar informações e dados

experimentais de sua pesquisa de doutorado.

Aos meus pais Marisa e Paulo pelo incentivo em buscar novos conhecimentos e

desafios e a Tatiana pelo apoio e paciência.

Resumo Alencar, Leandro de Souza; Naccache, Mônica Feijó. Simulação Numérica

de Reinício do Escoamento de Fluidos Tixotrópicos. Rio de Janeiro,

2016. 80p. Dissertação de Mestrado - Departamento de Engenharia

Mecânica, Pontifícia Universidade Católica do Rio de Janeiro.

A produção de petróleo em campos de água ultraprofunda requer linhas de

produção extensas no leito marinho onde a temperatura ambiente é baixa. O

conhecimento do comportamento tixotrópico do fluido produzido, aliado à

capacidade de simulação do problema de paradas de produção, pode permitir o

dimensionamento econômico de bombas e linhas para reiniciar a produção.

Diversos modelos que tentam explicar o comportamento tixotrópico de um fluido

existem na literatura. Os modelos que utilizam um parâmetro escalar para

representar a condição da estrutura do fluido são simples e possuem maior

aplicação prática, por isso, foram escolhidos. A equação de evolução para

tixotropia envolve uma equação diferencial não linear que deve ser acoplada à

resolução das demais equações de transporte necessárias. A utilização de uma

ferramenta de simulação numérica onde seja possível acoplar essas equações

permite ao engenheiro avaliar o problema para qualquer geometria e condição de

contorno. Neste trabalho é apresentado o algoritmo desenvolvido para permitir a

simulação tixotrópica no simulador comercial Fluent© para um fluido puramente

viscoso e incompressível. Dois cenários de repartida foram avaliados em detalhes

assim como o efeito da condição de contorno de entrada para a solução ao longo

do tempo do sistema. Uma relação do tempo para reinício com a pressão aplicada

é obtida para diferentes condições de tixotropia utilizando equações simples. A

simulação para o reinício do escoamento apresentou comportamento esperado e

qualitativamente semelhante aos resultados experimentais. Em conclusão, os

resultados apresentados comprovam o potencial uso da simulação numérica com

fluidos tixotrópicos especialmente para o problema de reinício permitindo o

dimensionamento de um sistema de produção de maneira eficiente.

Palavras-chave Tixotropia; simulação numérica; CFD; repartida de escoamento.

Abstract Alencar, Leandro de Souza; Naccache, Mônica Feijó. Numerical

Simulation of the Restart Problem for Thixotropic Fluids. Rio de

Janeiro, 2016. 80p. Dissertação de Mestrado - Departamento de Engenharia

Mecânica, Pontifícia Universidade Católica do Rio de Janeiro.

The offshore oil production in ultra-deep water fields requires long tie-backs

where seawater reaches low temperatures. The knowledge of the tixotropic

behavior of the produced fluid together with the simulation capacity to solve the

restart phenomenon, enables the economical design of pumps and lines with

reasonable precision and size to break the gel. Several mathematical models exist

in the literature that trys to explain this phenomenon, including scalar models in

which the structure of the fluid is represented by one single parameter. This type

of model is simple and with the most practical use. So it was chosen for this

thesis. A non-linear first order differential equation is used for the thixotropic

evolution relation and must be coupled with the standard transport equations. The

use of a numerical simulator that can solve this problem gives the engineer the

tool required to solve for any geometry the restart problem. In this thesis the

algorithm developed to implement the thixotropic simulation within the CFD

commercial tool Fluent© is presented for a purely viscous incompressible fluid.

Two restart cases are discussed in detail and the effect of the inlet boundary

condition is presented. A relationship between the restart time and the inlet

pressure is obtained for diferent thixotropic conditions. The simulation of the

restart problem showed an expected behaviour for thixotropic fluids with

qualitative similarities with experimental data. In conclusion, the results herein

presented prove the potencial use of CFD with thixotropic fluids enabling the

optimized design of production systems.

Keywords Thixotropy; numerical simulation; CFD; restart problem.

Sumário

1 Introdução 15

2 Revisão Bibliográfica 17 2.1 Fluido Não-Newtoniano e Tixotropia 17 2.1.1 Fluidos Puramente Viscosos e Independentes do Tempo 17 2.1.2 Fluido Viscoelástico 24 2.1.3 Fluidos Puramente Viscosos e Dependentes do Tempo 26 2.2 Trabalhos Anteriores 33

3 Formulação Matemática 36 3.1 Equações Fundamentais de Transporte 36 3.2 Discretização da Equação de Transporte Genérica 38 3.2.1 Funções de Interpolação para o Termo Advectivo 40 3.2.2 Acomplamento da Pressão e Velocidade 42 3.2.3 Discretização da Equação de Evolução da Tixotropia 42

4 Metodologia 45 4.1 Algoritmo para Resolução das Equações de Transporte no Fluent© 45 4.2 Principais Configurações no Fluent© 48 4.3 Teste de Independência de Malha 49

5 Resultados e Discussão 52 5.1 Propriedades Reológicas do Fluido Utilizado 52 5.2 Análises Qualitativas para o Reinício 53 5.2.1 Pressão de Entrada Baixa 53 5.2.2 Pressão de Entrada Alta 59 5.3 Avaliação das Configurações Iniciais 61 5.3.1 Condição de Contorno na Entrada para Lambda 61 5.3.2 Análise do Comprimento da Malha 64 5.4 Análises de Sensibilidade 68 5.4.1 Sensibilidade à Pressão de Entrada 68 5.4.2 Sensibilidade ao Parâmetro b 70 5.4.3 Sensibilidade ao Tempo de Equilíbrio 73 5.5 Validação Experimental 74

6 Conclusão 77

7 Referências bibliográficas 78

Lista de Figuras

Figura 1 Comportamento de diferentes tipos de fluidos, adaptado de Lima [3] 18

Figura 2 Viscosidade aparente e taxa de deformação para o modelo power-law (linha contínua) e viscosidade medida (linha pontilhada). 19

Figura 3 Comportamento dos Fluidos Viscoplásticos 21

Figura 4 Viscosidade por taxa de deformação para o modelo de Bingham e H-B 22

Figura 5 Viscosidade por tensão para o modelo de

Papanastasiou (linhas contínuas) para dois valores de η0 23

Figura 6 Viscosidade por tensão para o modelo utilizado no Fluent© (linha contínua) para dois valores de taxa crítica 23

Figura 7 Viscosidade por tensão para o modelo proposto por de Souza Mendes [10] 24

Figura 8 Análogo mecânico do modelo viscoelástico proposto por Maxwell e modelo viscoso 25

Figura 9 Comportamento viscoelástico – adaptado de Chhabra & Richardson [1] 25

Figura 10 Análogo mecânico do modelo de Jeffrey 26

Figura 11 Comportamento de viscosidade aparente dependente do tempo 27

Figura 12 Loop de histerese para diversos testes em sequência 28

Figura 13 Comportamento de viscosidade aparente dependente do tempo – adaptado deBarnes [12] 29

Figura 14 Teste com controle da tensão e o resultado na viscosidade 29

Figura 15 Análogo mecânico proposto por de Souza Mendes [9] 33

Figura 16 Resumo com as equações que devem ser resolvidas para o escoamento do fluido tixotrópico 38

Figura 17 Cálculo da variável nas faces 40

Figura 18 Interpolação no método de primeira ordem para obter

o valor da variável ϕl, referente a face entre as células A e B 40

Figura 19 Interpolação no método de primeira ordem, diferenças centrais e power-law para obter o valor da variável

ϕl, referente a face l entre as células A e B 41

Figura 20 Interpolação no método de segunda ordem e QUICK

para obter o valor da variável ϕl, referente a face entre as células A e B 42

Figura 21 Proposta de linearização para o termo fonte. Adaptado de [30] 43

Figura 22 Algoritmo para resolução das equações de transporte no Fluent©. Adaptado de [8] 45

Figura 23 Funções adicionais requeridas para implementar a tixotropia no Fluent© para o algoritmo segregado 46

Figura 24 Velocidade Axial média a 10cm da saída em função do número de células utilizadas na rede 50

Figura 25 Erros relativos à malha mais refinada (2160x180). 50

Figura 26 Perfil de velocidade obtida pelo cálculo analítico e no modelo do Fluent 51

Figura 27 Variação da estrutura λ do fluido para λeq igual a 0,5 em função dos parâmetros a e b 52

Figura 28 Cenário inicial apresentando a pressão aplicada na entrada e a área rachurada do duto em que a simulação é realizada (sistema axisimétrico) 53

Figura 29 Velocidade Média do Fluido 53

Figura 30 Mapa da estrutura do fluido para diferentes tempos ao longo do reinício com pressão de entrada de 1500 Pa e teq de 0,6 s 54

Figura 31 Perfil no raio para a estrutura do fluido em diferentes tempos com pressão de entrada de 1500 Pa e teq de 0,6 s 55

Figura 32 Taxa de deformação no centro em diferentes posições ao longo do comprimento 55

Figura 33 Estrutura do fluido no centro em diferentes instantes de tempo 56

Figura 34 Velocidade Média do Fluido para gradiente de λ igual a zero e λ igual a 0,2 como condição de contorno na entrada 57

Figura 35 Mapa da estrutura do fluido para diferentes tempos ao longo do reinício com condição de contorno na entrada para λ igual a 0,2 e teq de 0,6 s 58

Figura 36 Estrutura do fluido no centro em diferentes instantes de tempo para λ igual a 0,2 na entrada 59

Figura 37 Número de Reynolds e Comprimento de Entrada calculados próximo a saída ao longo do tempo para pressão aplicada alta (4500Pa) 59

Figura 38 Mapa da estrutura do fluido para diferentes tempos ao longo do reinício com pressão de entrada de 4500 Pa e teq de 0,6 s 60

Figura 39 Perfil de velocidade em diferentes instantes de tempo para pressão elevada aplicada 61

Figura 40 Velocidade média do fluido a 3,5 m da entrada para pressão de entrada de 1500 Pa e teq de 0,6 s para diferentes condições de entrada para λ 62

Figura 41 Posição média do fluido a partir da integral da velocidade média a 3,5m da entrada para diferentes condições de entrada para λ 63

Figura 42 Velocidade do fluido ao longo do comprimento para a pressão de entrada de 1500 Pa e teq de 0,6 s para diferentes condições de entrada para λ no tempo final 63

Figura 43 Estrutura do fluido para diferentes condições de entrada de λ no tempo final para a pressão aplicada de 1500 Pa e teq de 0,6 s 64

Figura 44 Velocidade média no tempo para as duas razões comprimento-diâmetro (velocidade medida na razão posição comprimento (x/L) igual a 0.975) 65

Figura 45 Velocidade Média em relação a velocidade média final para diferentes razões posição-comprimento para L/D 40 e L/D 20 65

Figura 46 Taxa de deformação no centro e velocidade média para pressão baixa aplicada e razão L/D igual a 40 e 20 66

Figura 47 Número de Reynolds e Comprimento de Entrada, para diferentes pressões aplicadas, calculados próximo a saída 67

Figura 48 Média do gradiente de velocidade axial ao longo do duto 67

Figura 49 Número de Reynolds medido na saída (x/L=0.975) para as maiores pressões aplicadas 68

Figura 50 Velocidade média do fluido em função do tempo e da pressão de entrada para teq igual a 0,6s 69

Figura 51 Velocidade média do fluido em função do tempo e da pressão de entrada para teq igual a 0,6s na escala logarítmica 69

Figura 52 Tempo para o reinício em função da pressão de entrada para dois tempos de equilíbrio 70

Figura 53 Derivada do termo tixotrópico (equação 40) para diferentes valores de b considerando λeq = 0.5 71

Figura 54 Derivada do termo tixotrópico (equação 40) para diferentes valores de a considerando λeq = 0.5 71

Figura 55 Velocidade média do fluido para a=0,5 e b igual a 0,5 e 2 72

Figura 56 Tempo para o reinício em função da pressão de entrada para dois valores de b 73

Figura 57 Velocidade média na saída (x/L=0.975) para diferentes tempos de equilíbrio 73

Figura 58 Tempo de repartida em função do tempo de equilíbrio para duas pressões de entrada (baixa e alta) 74

Figura 59 Pressão aplicada nos experimentos 75

Figura 60 Velocidade obtida a partir de dados experimentais 75

Figura 61 Velocidade média obtida pela simulação numérica 76

Figura 62 Comparativo entre dado experimental e simulado 76

Lista de Tabelas

Tabela 1 Algoritmo para o termo fonte 47

Tabela 2 Algoritmo para o termo de fluxo 47

Tabela 3 Algoritmo para o termo transiente 47

Tabela 4 Algoritmo para a viscosidade 47

Tabela 5 Algoritmo para manipular valores anteriores 47

Tabela 6 Fatores de relaxação utilizados 48

Tabela 7 Resíduos máximos (em valor absoluto) 49

Tabela 8 Diferentes malhas utilizadas 49

Tabela 9 Parâmetros Reológicos da Solução de Laponita Utilizada como Referência 52

Lista de Abreviação

m Coeficiente de consistência

Γ Coeficiente difusivo

L Comprimento do duto

k Constante de Bingham

𝛾 Deformação

D Diâmetro do duto

CFD Dinâmica de Fluidos Computacional

2-D Duas dimensões

𝜆 Estrutura do fluido

𝜆𝑒𝑞 Estrutura do fluido em equilíbrio (regime permanente)

𝑆𝜙 Fonte

H-B Herschel-Bulkley

n Índice de comportamento

𝜌 Massa específica

FVM Método dos volumes finitos

G Módulo de cisalhamento

Pe Número de Peclet

x Posição ao longo do comprimento

p Pressão aplicada

R Raio do duto

L/D Razão comprimento por diâmetro do duto

r Posição ao longo do raio do duto

�̇� Taxa de deformação

�̇�𝑐 Taxa de deformação crítica

𝑡𝑒𝑞 Tempo de equilíbrio

α Tempo de relaxamento

𝜋 Tensor tensão

𝜏 Tensor tensão viscosa

𝜏 Tensão

𝜏0 Tensão limite de escoamento

𝜙 Variável dependente

�⃗� Vetor velocidade

𝜂 Viscosidade aparente

𝜂𝑣 Viscosidade aparente (modelo de tixotropia)

𝜂𝑒𝑞 Viscosidade aparente em equilíbrio

𝜂0 Viscosidade aparente para cisalhamento nula

𝜂∞ Viscosidade aparente para cisalhamento muito alta

15

1 Introdução

O escoamento de fluidos tixotrópicos é um tema relevante em diversas áreas

da indústria. O transporte de tintas, alimentos em pasta, etc. por dutos rígidos e

linhas flexíveis dentro de uma instalação industrial requer o conhecimento do

comportamento tixotrópico destes fluidos para permitir dimensionar bombas,

tempos de parada máxima e outros de maneira econômica e eficiente. Na indústria

de petróleo a ocorrência de óleos gelificantes é também de extrema importância e

será a área alvo deste estudo.

A produção em campos de água ultraprofunda, cenário mais comum no

Brasil, implica em linhas de produção extensas no leito marinho que geram

desafios de garantia de escoamento diversos. Devido à baixa temperatura

ambiente no leito marinho, eventuais paradas de produção causam o resfriamento

dos fluidos dentro destas linhas. Especificamente para óleos parafínicos, com

tensão limite de escoamento elevada, este resfriamento promove a gelificação dos

fluidos. O reinício nestes casos pode ser possível somente com a aplicação de uma

elevada pressão pelo navio produtor para quebrar o gel formado. Desta forma, é

importante definir a pressão mínima para quebra do gel e o tempo necessário de

aplicação desta força. O conhecimento do comportamento tixotrópico deste

sistema aliado à capacidade de simulação do problema, pode permitir o

dimensionamento de bombas e linhas com classe de pressão compatível de

maneira segura e sem excessos.

Diversos modelos que tentam explicar o comportamento tixotrópico de um

fluido existem na literatura, desde sistemas complexos que descrevem o

comportamento das iterações intermolecures até sistemas algébricos simples que

modificam, por exemplo, a equação constitutiva de viscosidade utilizando um

parâmetro escalar para representar o nível estrutural do fluido. Estes modelos que

utilizam o parâmetro escalar são relativamente simples e possuem maior aplicação

prática. A dificuldade repousa em conseguir obter os parâmetros de ajuste

característicos para aquele fluido e/ou condição.

A equação de evolução para tixotropia a resolver é uma equação diferencial

de primeira ordem não linear e que deve ser acoplada à resolução das demais

equações de transporte necessárias para avaliar o escoamento de um fluido. A

16

utilização de uma ferramenta de simulação numérica onde seja possível acoplar

todas essas equações permite ao engenheiro avaliar, para qualquer geometria e

condição de contorno, um problema envolvendo um fluido tixotrópico. Desta

forma, é importante desenvolver uma metodologia para que esta ferramenta seja

capaz de utilizar uma equação de evolução, para representar corretamente o fluido

em questão.

O principal objetivo deste trabalho é desenvolver uma metodologia para

acoplar o fenômeno da tixotropia ao simulador numérico comercial Fluent© e

permitir a avaliação de escoamentos envolvendo fluidos tixotrópicos. Para isso, o

comportamento qualitativo do reinício de escoamento de um fluido qualquer será

avaliado e comparado com dados experimentais disponíveis. Além disso, é

proposta uma relação entre a pressão de entrada e alguns parâmetros reológicos

com o tempo necessário paro reinício.

Neste trabalho serão apresentados em uma breve revisão bibliográfica as

características de fluidos não newtonianos, em especial aqueles com

comportamento tixotrópico, e as equações matemáticas disponíveis que tentam

caracterizar esses fluidos. As equações de transporte necessárias e as técnicas de

discretização utilizadas para aplicar essas equações em uma ferramenta numérica

serão discutidas. Adicionalmente, a metodologia aplicada para permitir que um

simulador numérico comercialmente disponível seja capaz de representar o

comportamento tixotrópico será também descrita. Por fim, alguns resultados

experimentais serão apresentados com o objetivo de avaliar qualitativamente a

capacidade do modelo em representar o comportamento esperado.

17

2 Revisão Bibliográfica

2.1 Fluido Não-Newtoniano e Tixotropia

Em um fluido newtoniano a relação entre a tensão e a taxa de deformação é

linear e a viscosidade é normalmente representada por μ. A viscosidade depende

somente da pressão e da temperatura. Um fluido não newtoniano apresenta uma

relação não linear entre a tensão e a taxa de deformação. O termo η, denominado

de viscosidade aparente, relaciona a tensão e a taxa de deformação para fluidos

não newtonianos [1].

Quando a tensão 𝜏 é função somente da deformação 𝛾, o material é dito

puramente elástico (𝜏 = 𝜏(𝛾)). Por outro lado, quando depende somente da taxa

de deformação, é dito puramente viscoso (𝜏 = 𝜏(�̇�)). Já quando depende de

ambos, deformação e cisalhamento, é denominado viscoelástico (𝜏 = 𝜏(𝛾, �̇�)).

Outros materiais requerem relações mais complexas que dependem do tempo ou

da história de deformação do fluido. [2]

Desta forma, os fluidos não newtonianos podem ser agrupados em três

classes, conforme proposto em [1] e [2], a saber:

- Puramente viscosos e independentes do tempo;

- Viscoelásticos e dependentes do tempo;

- Puramente viscosos e dependentes do tempo.

Esta classificação, apesar de arbitrária, auxilia o estudo, a modelagem e o

entendimento dos fluidos. Um fluido real pode apresentar o comportamento de

dois ou dos três grupos, mas é geralmente possível enquadrá-lo em apenas um

destes grupos com base em sua característica dominante [1].

2.1.1 Fluidos Puramente Viscosos e Independentes do Tempo

Fluidos nos quais a taxa de deformação é determinada somente pela tensão

de cisalhamento instantânea, ou seja, que não dependem do histórico ou do tempo

de aplicação desta tensão, possuem viscosidade classificadas como independente

do tempo. A Figura 1 apresenta o comportamento da tensão pelo cisalhamento

para diferentes comportamentos: pseudoplástico, dilatante e viscoplástico.

18

Figura 1 Comportamento de diferentes tipos de fluidos, adaptado de Lima [3]

O comportamento pseudoplástico (ou shear-thinning) é caracterizado pela

redução da viscosidade aparente com o aumento da taxa de deformação e

apresenta comportamento newtoniano para taxas de cisalhamento muito elevadas

e muito baixas [1]. A viscosidade aparente tende para um valor 0, quando o

cisalhamento tende a zero, e tende para um valor para taxa muito alta,

conforme ilustrada na Figura 2.

Os fluidos dilatantes apresentam comportamento similar aos fluidos

pseudoplásticos, mas a viscosidade aumenta com o aumento do cisalhamento.

Conforme apresentado por Chhabra & Richardson, a fase líquida lubrifica as

partículas sólidas quando a taxa de deformação é baixa, o que reduz a viscosidade

aparente do material. Entretanto, o fluido se dilata com o aumento da taxa de

deformação e a fase líquida não é mais capaz de preencher todo o vazio entre as

partículas sólidas para, desta forma, evitar o contato direto entre as mesmas (reduz

o efeito da lubrificação). Este fenômeno resulta, portanto, em aumento de fricção

interna e aumento da tensão cisalhante e da viscosidade [1].

Fluidos com comportamento viscoplástico são aqueles que não se deformam

quando submetidos a tensões menores do que a tensão limite de escoamento (𝜏0).

O fluido pode apresentar comportamento elástico ou mover-se como corpo rígido.

Acima da tensão limite, a estrutura do fluido se quebra e o mesmo se comporta

como um fluido viscoso. A partir deste momento, a relação entre tensão e taxa de

deformação pode ser linear, como um fluido newtoniano, ou não linear, como um

fluido pseudoplástico [3].

𝜏

𝛾 ̇

Pseudoplástico

Newtoniano

Dilatante

Viscoplásticos

19

Diversos modelos matemáticos estão disponíveis na literatura para descrever

o comportamento dos fluidos não newtonianos puramente viscosos. A seguir são

apresentados alguns destes modelos.

2.1.1.1 Modelo Power-Law

No modelo power-law ou Ostwald e de Waele, proposto em 1925 [4], a

relação entre tensão 𝜏 e taxa de deformação �̇� é dada pela relação

𝜏 = 𝑚(�̇�)𝑛 1

onde o coeficiente de consistência (m) e o índice de comportamento (n) são

obtidos empiricamente.

Observa-se que o fluido assume comportamento dilatante quando n é maior

que a unidade, comportamento pseudoplástico quando n é menor que a unidade e

comportamento newtoniano quando n é igual a um. Desta forma, o índice de

comportamento n indica a semelhança do fluido a um comportamento

newtoniano. Já o parâmetro m está relacionado à viscosidade medida e indica a

resistência do fluido em escoar.

A Figura 2 compara a viscosidade aparente pela taxa de deformação de um

fluido qualquer (linha pontilhada) e do modelo power-law ajustado (linha

contínua).

Figura 2 Viscosidade aparente e taxa de deformação para o modelo power-law (linha

contínua) e viscosidade medida (linha pontilhada).

Adaptado de Chhabra & Richardson [1]

0,0001

0,0010

0,0100

0,1000

1,0000

10,0000

0,01 0,10 1,00 10,00 100,00 1000,00 10000,00

Vis

cosi

dad

e A

par

en

te (

Pa.

s)

Taxa de Cisalhamento (1/s)

𝜇 0

𝜇 ∞

20

O modelo power-law apresenta algumas limitações, como, por exemplo, não

ser capaz de representar corretamente a viscosidade para taxas de cisalhamento

muito baixas ou muito altas, os valores não tendem para 𝜇0 ou 𝜇∞ [4], conforme

ilustra a Figura 2. Uma grande vantagem deste modelo é permitir com facilidade

sua utilização em análises analíticas [2].

2.1.1.2 Modelo Carreau-Yasuda

O modelo Carreau-Yasuda [4], ao contrário do modelo power-law,

considera os desvios na viscosidade aparente para os extremos de taxa de

deformação conforme apresentado na Figura 2. O modelo foi inicialmente

proposto por Carreau em 1968 e modificado por Yasuda em 1981 conforme

equação 2.

𝜂 − 𝜂∞

𝜂0 − 𝜂∞= [1 + (𝜆�̇�)𝑎]

𝑛−1𝑎 2

O parâmetro 𝜆 do modelo é um índice de consistência e a é um valor

adimensional que ajusta a transição entre a região de taxa de deformação próxima

a zero (𝜂0) e a região logarítmica linear do modelo de power-law [4].

2.1.1.3 Modelo de Bingham

Os modelos de power-law ou Carreau-Yasuda não descrevem o

comportamento viscoplástico. Já o modelo de Bingham, proposto em 1922 [2],

representa o fluido como sólido para tensões abaixo da tensão limite de

escoamento e como puramente viscoso acima desta tensão. Esse comportamento é

caracterizado pela equação 3.

𝜂(�̇�) = {

𝜇 +𝜏0

�̇� 𝑝𝑎𝑟𝑎 𝜏 ≥ 𝜏0

∞ 𝑝𝑎𝑟𝑎 𝜏 < 𝜏0

3

O modelo de Bingham é comumente utilizado para representar fluidos de

perfuração na indústria de petróleo, onde a viscosidade 𝜇 e a tensão limite 𝜏0 são

funções da fração de partículas sólidas no fluido, do diâmetro destas partículas e

da viscosidade da fase líquida [2].

Este modelo é capaz de representar o comportamento linear entre tensão e

taxa de deformação de um fluido viscoplástico a partir da tensão limite de

escoamento, mas não é possível modelar fluidos com comportamento dilatante ou

21

pseudoplástico. A Figura 3 ilustra o modelo de Bingham e sua caracterísitica

newtoniana para tensão maior que a tensão limite.

Figura 3 Comportamento dos Fluidos Viscoplásticos

2.1.1.4 Modelo Herschel-Bulkley

Em 1926, Herschel e Bulkley propuseram uma generalização do modelo de

Bingham, seguindo a mesma lógica do modelo power-law de 1925 [5]. As três

curvas da Figura 3 são possíveis de representar pelo modelo Herschel-Bulkley (H-

B). A relação entre tensão 𝜏 e taxa de deformação é dada pela relação

𝜏 = {

𝜏0 + 𝑚�̇�𝑛 0

𝑝𝑎𝑟𝑎 𝜏 ≥ 𝜏0 𝑝𝑎𝑟𝑎 𝜏 < 𝜏0

4

onde os parâmetros m e n são os mesmos do modelo power-law (item 2.1.1.1 ).

Quando n=1 a equação 4 equivale à equação 3 (modelo de Bingham).

A viscosidade aparente é obtida partir da relação 𝜏 = �̇�η, para tensão maior

que a tensão limite de escoamento (𝜏 ≥ 𝜏0), conforme equação 5. A Figura 4

representa a viscosidade pela taxa de deformação para um fluido representado

pelo modelo de Bingham e pelo modelo de Herschel-Bulkley (H-B).

𝜂 =𝜏0

�̇�+ 𝑚�̇�𝑛−1 𝑝𝑎𝑟𝑎 𝜏 ≥ 𝜏0 5

𝜏

𝛾 ̇

Bingham

𝜏0

H-B (n>1)

H-B (n<1)

22

Figura 4 Viscosidade por taxa de deformação para o modelo de Bingham e H-B

2.1.1.5 Modelos Regularizados de Viscosidade

Conforme a Figura 4, a viscosidade no modelo de Herschel-Bulkley (H-B)

tende a infinito para taxas de cisalhamento que tendem a zero. Em uma simulação

numérica, este comportamento deve ser aproximado por um modelo que limite o

valor máximo da viscosidade. Além disso, o comportamento de sólido, para

tensões abaixo da tensão limite de escoamento, é, na verdade, muitas vezes

representado pelo comportamento viscoso com viscosidade muito alta,

representada geralmente por 𝜂0. Desta forma, a equação de viscosidade deve ser

suavizada.

Diversas equações que tratam essa questão estão disponíveis na literatura.

Papanastasiou propôs em seu artigo [6], a partir da análise de dados empíricos de

fluidos com comportamento de Bingham quase ideal, a equação 6. Verifica-se que

esta equação não é singular, mesmo para taxa de deformação nula, ao aplicar a

regra de L’Hôpital.

𝜂 = 𝑚 +𝜏0

�̇�(1 − 𝑒−𝜂0�̇�) 6

A equação proposta por Papanastasiou pode ser modificada para adequar-se

ao comportamento de fluidos do tipo H-B conforme proposto em [7] e

apresentada na equação 7.

𝜂 = 𝑚�̇�𝑛−1 +𝜏0

�̇�(1 − 𝑒−𝜂0�̇�) 7

A Figura 5 apresenta a equação de H-B original (equação 5) e a formulação

modificada de Papanastasiou (equação 7).

1E-06

1E-04

1E-02

1E+00

1E+02

1E+04

1E+06

1E+08

1E-06 1E-03 1E+00 1E+03 1E+06 1E+09

Vis

cosi

dad

e A

par

en

te (

Pa.

s)

Taxa de cisalhamento (1/s)

Bingham

Herschel-Bulkley

τ0 = 5Pa

mH-B ou µbingham = 2 Pa.s nH-B = 0.5

23

Figura 5 Viscosidade por tensão para o modelo de Papanastasiou (linhas contínuas) para

dois valores de 𝛈𝟎

O simulador computacional de dinâmica de fluidos (CFD) Fluent© possui

outra proposta de regularização dada pelas relações a seguir [8], onde a taxa

crítica �̇�𝑐 deve ser informada pelo usuário:

Se �̇� > �̇�𝑐:

Se �̇� ≤ �̇�𝑐:

𝜂 =𝜏0

�̇�+ 𝑚�̇�n−1

𝜂 =𝜏0(2 − �̇�/�̇�𝑐 )

�̇�𝑐+ 𝑚�̇�𝑐

n−1[(2 − 𝑛) + (𝑛 − 1)�̇�/�̇�𝑐]

8

A taxa crítica possui o mesmo efeito que a viscosidade 𝜂0 no modelo de

Papanastasiou, comparando a Figura 6 à Figura 5. Por existir duas equações, deve-

se garantir que a taxa de deformação crítica é baixa o suficiente para evitar

descontinuidades grandes ao trocar de equação.

Figura 6 Viscosidade por tensão para o modelo utilizado no Fluent© (linha contínua) para

dois valores de taxa crítica

1,E-04

1,E-02

1,E+00

1,E+02

1,E+04

1,E+06

1,E+08

1,E-02 1,E-01 1,E+00 1,E+01 1,E+02 1,E+03 1,E+04

Vis

cosi

dad

e a

par

en

te (

Pa.

s)

Tensão (Pa)

Herschel-Bulkley

η0 = 1e6

η0 = 1e5

1,E-05

1,E-03

1,E-01

1,E+01

1,E+03

1,E+05

1,E+07

1,E-03 1,E-01 1,E+01 1,E+03 1,E+05 1,E+07

Vis

cosi

dad

e a

par

en

te (

Pa.

s)

Tensão (Pa)

Herschel-Bulkley

γc = 1e-5

γc = 1e-6

24

Em 2004, de Souza Mendes e Dutra [7] apresentaram uma proposta de

viscosidade semelhante à de Papanastasiou, dada pela equação 9. A viscosidade é

igual a 𝜂0 para taxa de deformação nula enquanto que, para Papanastasiou,

tenderia a um valor infinito.

𝜂 = 𝑚�̇�𝑛−1 +

𝜏0

�̇�(1 − 𝑒

−𝜂0�̇�𝜏0

⁄ ) 9

A equação anterior foi modificada em trabalhos posteriores por de Souza

Mendes [9] e outros [10] para a equação 10. A equação para viscosidade proposta

tende para 𝜂0 quando a taxa de deformação é nula, assim como na equação 9, e

tende para η∞ quando a taxa tende para infinito (vide Figura 7). Essa

característica de ser limitada nos dois extremos é importante para evitar erro de

ponto flutuante em simulação numérica conforme discutido em [11].

𝜂 = (𝑚�̇�𝑛−1 +

𝜏0

�̇�) (1 − 𝑒

−𝜂0�̇�𝜏0

⁄ ) + 𝜂∞ 10

Figura 7 Viscosidade por tensão para o modelo proposto por de Souza Mendes [10]

2.1.2 Fluido Viscoelástico

Os fluidos apresentados na seção anterior (item 2.1.1) possuem

comportamento reológico puramente viscoso, onde a viscosidade é dependente da

taxa de deformação. Por outro lado, alguns fluidos apresentam um comportamento

elástico pronunciado [1]. Estes fluidos, classificados como viscoelásticos, são

capazes de absorver e recuperar parte da tensão cisalhamento aplicada.

1,E-05

1,E-03

1,E-01

1,E+01

1,E+03

1,E+05

1,E+07

1,E-03 1,E-01 1,E+01 1,E+03 1,E+05 1,E+07

Vis

cosi

dad

e a

par

en

te (

Pa.

s)

Tensão (Pa)

Herschel-Bulkley

de Souza Mendes

25

O modelo proposto por Maxwell [4] foi um dos primeiros a tentar explicar

este padrão ao unir o comportamento de sólidos elásticos (lei de Hooke) ao

comportamento de fluidos viscosos, a partir da equação

𝜏 +

𝜇

𝐺

𝜕𝜏

𝜕𝑡= −𝜇�̇� 11

onde G é o módulo de cisalhamento, e o termo 𝜇/𝐺 é o de tempo de relaxamento

representado por α [4].

Em regime permanente o modelo de Maxwell indica que a tensão é

proporcional a taxa de deformação, como um fluido Newtoniano. Quanto maior a

taxa de variação no tempo da tensão, maior será a contribuição do efeito elástico

apresentado na equação 11.

A Figura 8 compara o análogo mecânico do modelo proposto por Maxwell

ao modelo puramente viscoso apresentado no item 2.1.1.

Figura 8 Análogo mecânico do modelo viscoelástico proposto por Maxwell e modelo viscoso

O experimento proposto por Chhabra & Richardson [1] ilustra o

comportamento de fluido viscoelásticos e o efeito do tempo de relaxamento λ: a

partir do fluido em repouso, aplica-se um deslocamento instantâneo e constante

(taxa de deformação nula) e mede-se a tensão. Os resultados são apresentados na

Figura 9. Quanto maior o tempo de relaxamento, mais o fluido viscoelástico se

aproxima da resposta puramente elástica.

Figura 9 Comportamento viscoelástico – adaptado de Chhabra & Richardson [1]

τ

μ G

η

τ Viscoelástico de Maxwell

Viscoso

log𝜏

𝑡(𝑠)

λ=0.5

λ=1

λ=5

Aumento do

comportamento

elástico

Aumento do

comportamento

viscoso

26

Os fluidos viscoelásticos possuem efeito de memória, onde a tensão em um

dado instante t depende da taxa de deformação aplicada no mesmo instante e das

taxas dos instantes anteriores. A partir da integral da equação 11 é possível obter a

equação 12 [4], onde é claro esse fenômeno de memória. O efeito das taxas em

tempos anteriores t’ é cada vez menor quanto mais afastado do instante atual t,

conforme representado pelo termo exponencial.

𝜏 = −∫ (𝐺𝑒−(𝑡−𝑡′)

𝜇𝐺) �̇�(𝑡′)𝑑𝑡′

𝑡

−∞

12

Outros modelos podem ser propostos ou estão disponíveis na literatura [1, 2,

4, 5]. A proposta de Jeffrey [4], por exemplo, considera também a aceleração

(segunda derivada da deformação 𝛾) conforme equação 13 e análogo mecânico

representado na Figura 10.

𝜏 +

𝜂2

𝐺

𝜕𝜏

𝜕𝑡= −𝜇 (�̇� +

𝜂2

𝐺

𝜕�̇�

𝜕𝑡) 13

Figura 10 Análogo mecânico do modelo de Jeffrey

2.1.3 Fluidos Puramente Viscosos e Dependentes do Tempo

No item 2.1.1 foram apresentados os fluidos cuja tensão depende do produto

da taxa de deformação pela viscosidade aparente. Esta viscosidade pode ser

constante (fluido Newtoniano) ou variável, a depender da taxa de por exemplo. Já

no item 2.1.2, o efeito elástico é considerado e uma equação diferencial de

primeira ordem é apresentada para relacionar a tensão com a taxa de deformação e

a viscosidade aparente. Esta viscosidade, assim como nos fluidos puramente

viscosos (item 2.1.1), não varia com o tempo. Outro tipo de comportamento não

Newtoniano é denominado tixotropia. Neste caso, a viscosidade do fluido depende

do tempo de cisalhamento, e o comportamento é reversível. O inverso deste

comportamento também pode ocorrer, i.e., a viscosidade aumentando com o

tempo de cisalhamento, e, neste caso, o fluido é denominado anti-tixotrópico.

τ

η1

η2

G

27

Fluidos tixotrópicos tem como característica a quebra da sua estrutura ao

longo do tempo com a aplicação de cisalhamento e, em repouso, essa estrutura se

regenera lentamente. A quebra dessa estrutura implica em redução da viscosidade

aparente. Estes fluidos são capazes de retornar ao seu estado inicial após a

aplicação de uma taxa ou tensão de cisalhamento. O cisalhamento pode também

resultar em um comportamento inverso ao do fluido tixotrópico, ou seja, causa

crescimento da estrutura ao aumentar a colisão entre as partículas que, quando se

aproximam, se atraem e aderem-se uma as outras. Este comportamento é dito anti-

tixotrópico [12]. A Figura 11 ilustra os dois comportamentos.

Figura 11 Comportamento de viscosidade aparente dependente do tempo

O termo tixotropia foi sugerido por Peterfi em 1927 ao unir duas palavras

gregas: thixis (agitar) e trepo (mudar) [13]. Desde o início do século 20 o

comportamento tixotrópico é discutido, mas com significativo avanço a partir de

1950. Os primeiros trabalhos nessa área iniciaram com Schalek e Szegvari que,

em 1923, verificaram que géis aquosos com Fe2O3 poderiam ser transformados

em soluções líquidas quando agitados e, após um tempo em repouso, essas

soluções retornavam ao estado gelificado [12]. A definição de tixotropia evoluiu

desde então: de mudanças reversíveis nas propriedades de um fluido escoando

para um estado gelificado quando em repouso, para incluir a dependência do

tempo com a duração da aplicação de uma tensão cisalhante em 1967 [14].

A definição atual, conforme a IUPAC (International Union of Pure and

Applied Chemistry) é: a redução da viscosidade de um fluido inicialmente em

repouso com a aplicação de um cisalhamento ao longo do tempo. Se a viscosidade

é recuperada com o tempo após a aplicação deste cisalhamento, o fluido é dito

tixotrópico [15].

τ

�̇�

Tixotropia

Anti-tixotropia

28

Com essa definição, é introduzida a dependência do tempo, ou efeito de

memória. Além disso, é importante destacar que o fluido tixotrópico deve ser

reversível e apresenta geralmente histerese conforme será visto no item 2.1.3.1.

2.1.3.1 Caracterização da Tixotropia

A caracterização da tixotropia de um fluido é complexa, pois envolve

dependência do tempo e do histórico, mudanças irreversíveis, histerese e

sedimentação [12, 14]. A seguir são apresentadas algumas técnicas para avaliar as

propriedades da tixotropia em uma dada amostra.

A técnica do loop de histerese foi proposta por Green e Wetlmann em 1943

[16] e consiste em aplicar uma taxa de deformação crescente por um tempo e

reduzir gradativamente essa taxa. Desta forma, é possível obter o gráfico de

tensão por taxa de deformação, apresentado ilustrativamente na Figura 12, para

diversos ciclos realizados em sequência.

A forma da curva deste gráfico varia a depender do tipo de fluido, tempo de

repouso até o início do experimento, etc. [17]. Embora este teste seja simples, a

taxa de deformação e o tempo estão acoplados e, desta forma, não é possível

separar cada um destes efeitos no comportamento do fluido [12].

Figura 12 Loop de histerese para diversos testes em sequência

A variação em degraus da taxa de deformação até atingir o regime

permanente, permite verificar o efeito do cisalhamento no comportamento do

fluido e separar o efeito do tempo. A reversibilidade do fluido é possível verificar

executando os passos na ordem inversa [17]. O aumento e redução da viscosidade

ao longo do tempo com a variação da taxa é mostrado na Figura 13, que

exemplifica o teste para degraus de taxa de deformação.

τ

�̇�

τ

�̇�

τ

�̇�

29

Figura 13 Comportamento de viscosidade aparente dependente do tempo – adaptado

deBarnes [12]

O controle da tensão, ao invés da taxa de deformação, permite avaliar o

comportamento do fluido quando esta tensão está abaixo da tensão limite, por

exemplo. A redução em degrau da tensão aplicada pode resultar em um leve

aumento da viscosidade até outro valor em regime permanente ou o aumento até o

fluido parar de escoar. Os comportamentos da tensão e viscosidade são ilustrados

na Figura 14. Este teste é interessante para validar modelos de tixotropia em

ensaios laboratoriais [17].

Figura 14 Teste com controle da tensão e o resultado na viscosidade

Adaptado de Mewis & Wagner [17]

Outro teste realizado é aquele em que o fluido, inicialmente em repouso, é

sujeito a uma tensão ou taxa cisalhante constante, como em um cenário de

repartida do escoamento. De acordo com Barnes [12] a maioria dos fluidos

tixotrópicos quando em repouso, apresentam comportamento viscoelástico. Desta

forma, é esperado que no teste de repartida, para um degrau na taxa de

deformação, a tensão suba rapidamente atingindo um pico de máximo e, em

seguida, reduza gradativamente (conforme Figura 13). O gráfico desta tensão de

pico por tempo de repouso pode indicar o tempo de recuperação do fluido após o

cisalhamento [17].

τ

tempo

�̇�

tempo

breakdown

buildup

Permanente

t1 t2 t1 t2

η

tempo

τ

tempo

1

τ0

2

1

2

30

2.1.3.2 Mecanismos da Tixotropia

A tixotropia pode ser explicada a partir da microestrutura do fluido, criada

como resultado das forças de atração entre as partículas deste fluido. A

microestrutura pode estar associada, por exemplo, a floculação de partículas,

alinhamento de fibras, distribuição espacial das partículas ou a associação

molecular em soluções poliméricas [12].

Essas forças de união são capazes de gerar uma rede de partículas que

podem ser rompidas pelas forças hidrodinâmicas. Desta forma, quando inicia o

escoamento, ou a aplicação de uma tensão, a estrutura se quebra em estruturas

menores ou flóculos [17]. O aumento da taxa de deformação implica em estruturas

cada vez menores. Quando a tensão é reduzida ou encerrada, as forças de adesão

permitem que a rede de partículas se reconstrua.

Quando há equilíbrio entre a taxa de quebra da estrutura, causada pelas

forças de cisalhamento, e entre a taxa de formação de novas conexões, causada

pela união dos flóculos quando se colidem, o sistema está em regime permanente

[18].

Já o efeito de histerese (conforme apresentado no item 2.1.3.1) pode ocorrer

porque diferentes subestruturas, ou flóculos, podem ser formadas a depender da

taxa e da intensidade do cisalhamento sofrido. Essas subestruturas podem se unir a

outras formando subestruturas maiores com a redução do cisalhamento e

dependem novamente da taxa e da intensidade desta redução [17]. Observa-se,

portanto, que o caminho percorrido na quebra da microestrutura do fluido pode ser

diferente do caminho de formação, o que explica a histerese na viscosidade do

fluido, por exemplo.

2.1.3.3 Modelos para a Tixotropia

Existem na literatura diversos modelos matemáticos que representam, de

alguma forma, os fenômenos característicos da tixotropia. Estes modelos podem

ser divididos em três grandes grupos, conforme proposto em [12, 14], a depender

da estratégia adotada para descrever o comportamento da microestrutura do

fluido, a saber:

31

i. Modelos fenomenológicos, que utilizam propriedades reológicas que

dependem de funções de memória;

ii. Modelos cinéticos, que descrevem de maneira direta o comportamento da

microestrutura;

iii. Modelos indiretos, que utilizam um parâmetro escalar para descrever a

microestrutura.

A estratégia básica adotada pelos modelos fenomenológicos consiste em

definir a viscosidade ou tensão como uma função que dependa diretamente do

tempo atual e do passado [14, 18]. Estes modelos foram primeiramente propostos

por Slibar & Paslay [18] a partir da equação de Bingham adaptada com a tensão

limite de escoamento variável em função da taxa de deformação e tempo. Estes

modelos descrevem o comportamento do fluido sem relacionar diretamente o

estado da microestrutura do mesmo.

Com o objetivo de descrever diretamente a microestrutura do fluido, os

modelos cinéticos foram propostos [17]. Goldberg [18] propôs descrever o

número de conexões entre as partículas do fluido pela equação diferencial

−

𝑑𝑛

𝑑𝑡= 𝑘1𝑛�̇�𝑎 − 𝑘2(𝑛0 − 𝑛)�̇�𝑏 14

onde 𝑛0 indica a quantidade de conexões no fluido totalmente estruturado, 𝑛 o

número de conexões ativas, 𝑎 e 𝑏 correspondem respectivamente ao expoente

para a taxa de quebra e construção da estrutura.

A partir do número de conexões ativas, pode-se, por exemplo, calcular o

volume efetivo ocupado pelos flóculos e assim o impacto na viscosidade com base

no arraste hidrodinâmico causado por este volume, conforme descrito por Mewis

[17] e proposto por Potanin e cols. [19].

O terceiro grupo de modelos são aqueles que utilizam um parâmetro escalar

que define a estrutura para um dado instante de tempo. Esta estratégia se mostra

mais simples [14, 18], e será a adotada nos estudos aqui apresentados. O conceito

deste parâmetro, representado pelo símbolo λ, foi inicialmente proposto por

Moore em 1959 [20] e segue a relação apresentada na equação 15, com a e b

iguais a um. O parâmetro λ representa a fração de conexões ativas e guarda,

portanto, relação direta com 𝑛 do modelo microestrutural proposto por Goodeve

(equação 14).

32

𝑑𝜆

𝑑𝑡= 𝑘1𝜆�̇�

𝑎 − 𝑘2(1 − 𝜆)�̇�𝑏 15

Os modelos cinéticos e indiretos precisam relacionar a microestrutura do

fluido, descrita pelas equações diretas ou indiretas propostas, com seu

comportamento reológico [17]. Portanto, uma relação entre o parâmetro de

estrutura (equação 15) ou o número de conexões (equação 14) com uma ou mais

propriedades reológicas deve ser proposta. O modelo de Moore, e a maioria dos

modelos propostos na literatura, relaciona a viscosidade aparente com o parâmetro

estrutural [18]. A viscosidade proposta por Moore é representada por

𝜂 = 𝜂∞ + (𝜂0 − 𝜂∞)𝜆 16

onde 𝜂0 e 𝜂∞ são os valores extremos que a viscosidade pode assumir

para, respectivamente, valores muito baixos e muito altos de taxa de deformação

[21].

Com o objetivo de representar fluidos com tensão limite de escoamento,

Hooska [22] associou o parâmetro construtivo 𝜆 na equação de Herschel-Bulkley

(vide item 2.1.1.4 ). A tensão limite de escoamento passou a ser representada por

uma tensão 𝜏0 constante e outra tensão 𝜏0𝜆 devido à tixotropia, , conforme

equação 17, obtida de [23].

𝜏 = {

(𝜏0 + 𝜆𝜏0𝜆) + (𝑘 + 𝜆 𝛥𝑘) �̇�𝑛 0

𝑝𝑎𝑟𝑎 𝜏 ≥ 𝜏0 + 𝜆𝜏0𝜆 𝑝𝑎𝑟𝑎 𝜏 < 𝜏0 + 𝜆𝜏0𝜆

17

onde k corresponde a viscosidade permanente e Δ𝑘 a viscosidade em função da

tixotropia.

O parâmetro 𝜆 no modelo proposto por Hooska é dado pela equação 18,

onde 𝜆 varia de 1 (fluido totalmente estruturado) até 0 (totalmente sem estrutura).

A semelhança da proposta de Hooska (equação 18) com a generalização do

modelo inicialmente proposto por Moore (equação 15) é clara. Oito parâmetros

são necessários para definir as relações propostas, a saber: 𝜏0, 𝜏0𝜆, 𝑘, Δ𝑘, 𝑛, 𝑘1, 𝑎

e 𝑘2. No trabalho de Lima [3] são descritos os experimentos reológicos

necessários para a caracterização destes parâmetros.

𝑑𝜆

𝑑𝑡= 𝑘1𝜆�̇�

𝑎 − 𝑘2(1 − 𝜆) 18

de Souza Mendes propôs um modelo robusto [24] que inicialmente

considera também o comportamento elástico do fluido com base no modelo

viscoelástico proposto por Maxwell (descrito no item 2.1.2 e Figura 8) e em 2011

[25] atualizou a relação proposta para o modelo de Jeffrey (vide Figura 10). O

33

análogo mecânico proposto por Jeffrey e modificado com a constante elástica e

viscosa como funções do parâmetro construtivo 𝜆 é apresentado na Figura 15.

Figura 15 Análogo mecânico proposto por de Souza Mendes [9]

A equação constitutiva para a tensão é, desta forma, obtida com base nos

deslocamentos e tensões consideradas na Figura 15 e segue a equação 19, que é

semelhante a obtida para o modelo viscoelástico de Jeffrey.

𝜏 = 𝜂𝑣(�̇� + 𝜃2�̈�) − 𝜃1�̇� 19

onde 𝜃1 =𝜂𝑠

𝐺𝑠 , 𝜃2 = 𝜃1

𝜂∞

𝜂𝑠+𝜂∞ e 𝜂𝑣 = (𝜂𝑠 + 𝜂∞).

Considerando que a viscosidade aparente 𝜂𝑣 deve variar entre um valor

muito alto 𝜂0, fluido estruturado (𝜆 = 1), e muito baixo 𝜂∞ (𝜆 = 0) e que o

módulo de elasticidade deve aumentar com a redução da estrutura do fluido [9]

(redução de 𝜆), as relações da equação 20 são propostas para a viscosidade

aparente e o módulo de elasticidade.

𝜂𝑣(𝜆) = (

𝜂0

𝜂∞)𝜆

𝜂∞

𝐺𝑆(𝜆) =𝐺𝑆

𝜆𝑚

20

O parâmetro de construção 𝜆 pode ser obtido a partir da equação de

evolução 21 proposta por de Souza Mendes [26].

𝑑𝜆

𝑑𝑡=

1

𝑡𝑒𝑞[(1 − 𝜆)𝑎 − (1 − 𝜆𝑒𝑞)

𝑎(

𝜆

𝜆𝑒𝑞)

𝑏

] 21

onde o tempo de equilíbrio 𝑡𝑒𝑞, 𝑎 e 𝑏 são constantes características do fluido e 𝜆𝑒𝑞

corresponde ao parâmetro construtivo em equilíbrio (regime permanente).

2.2 Trabalhos Anteriores

A modelagem numérica do comportamento de fluidos não newtonianos, em

especial aqueles com comportamento tixotrópico, utilizando volumes finitos

τ

η∞

Gs(λ)

γ

γe γv

η𝑠(λ)

34

(FVM, do inglês finite-volume method) foi realizada em diversos trabalhos.

Alguns destes trabalhos são apresentados a seguir.

O deslocamento de um fluido tixotrópico em um duto por um fluido

newtoniano foi avaliado por de Souza Mendes e colaboradores [26]. O algoritmo

numérico proposto é capaz de avaliar a posição da interface entre os fluidos ao

longo do tempo, sem considerar compressibilidade ou elasticidade. A interface

entre os fluidos é assumida plana e se move conforme a velocidade média do

escoamento. A equação constitutiva utilizada para a tixotropia segue a equação

21, com a e b iguais a um, o que permite obter uma relação linear entre o

parâmetro tixotrópico e a taxa de deformação. O comportamento do escoamento é

comparado para diferentes parâmetros, por exemplo, tempo de equilíbrio e

viscosidade do fluido newtoniano entrante. A proposta numérica apresentada

também permitiu verificar o comportamento de avalanche.

Wachs e colaboradores [23] propõem um modelo FVM 1.5D capaz de

representar o comportamento tixotrópico para fluido viscoplástico compressível

em problemas de repartida. No modelo 1.5-D, apenas a componente axial da

velocidade é considerada, apesar dela depender tanto da posição axial como

radial. O modelo para a tixotropia de Houska foi utilizado, mas simplificado para

obter uma relação linear entre variação de λ e taxa de deformação. Outra

simplificação realizada foi considerar que apenas a tensão limite de escoamento

varia com o parâmetro construtivo (a viscosidade é mantida constante). Wachs e

cols. mostram que, ao considerar a tixotropia e a compressibilidade de um fluido

viscoplástico em um tubo, é possível, em alguns casos, reiniciar o escoamento

para pressões menores que 4𝜏0𝐿/𝐷, onde 𝜏0 é a tensão limite, L o comprimento e

D o diâmetro do tubo. [23]

No trabalho de Ahmadpour & Sadeghy [24] é utilizado o método de FVM

para modelar o reinício da produção de um óleo gelificante, conforme algoritmo

desenvolvido por Ahmadpour em sua tese [13]. Com base no modelo 2-D de um

tubo, diferentes análises de sensibilidade são realizadas para um fluido que segue

as equações constitutivas de viscoelasticidade propostas por Dullaert & Mewis

[27]. Diferentes análises de sensibilidade são realizadas, variando parâmetros de

elasticidade, estrutura e compressibilidade, para verificar a influência destes

parâmetros no tempo de repartida até o regime permanente.

35

O trabalho de Potanin [28] apresenta uma transição para a utilização de

ferramentas computacionais disponíveis comercialmente. Um reômetro é

representado por uma ferramenta de CFD com o objetivo de reproduzir

experimentos reológicos com fluidos tixotrópicos. Uma relação linear entre a

tixotropia e a taxa de deformação é considerada, para simplificar a solução do

problema. Os resultados apresentados mostram boa coerência com os

experimentos de caracterização realizados e demonstra a validade em utilizar

ferramentas comerciais de CFD, tais como Fluent© ou Flow-3D©, na modelagem

de fluidos tixotrópicos. [28]

Além do trabalho de Potanin [28], não foram encontrados outros trabalhos

disponíveis na literatura que utilizam estas ferramentas comerciais para

escoamento de fluido tixotrópico.

Desta forma, a proposta do presente estudo, desenvolver uma metodologia e

simular o comportamento tixotrópico de um fluido em uma ferramenta de CFD,

ainda não foi tema disponível na literatura e poderá servir de insumo para

trabalhos futuros.

36

3 Formulação Matemática

3.1 Equações Fundamentais de Transporte

O escoamento de um fluido é regido pela equação de conservação de massa,

de quantidade de movimento e de energia [4]. O problema a ser modelado será o

escoamento isotérmico do fluido com o objetivo de simplificar a análise e já que é

uma hipótese realista considerar que o fluido em repouso estaria em equilíbrio

térmico com o meio ambiente. Desta maneira, somente as equações de

conservação de massa e de quantidade de movimento são necessárias para a

obtenção da solução. As equações de conservação de massa e de quantidade de

movimento são apresentadas respectivamente pelas equações 22 e 23.

𝜕𝜌

𝜕𝑡= −∇. 𝜌�⃗� 22

𝜕

𝜕𝑡(𝜌�⃗� ) = −𝜌�⃗� . ∇�⃗� − ∇⃗⃗ . 𝜋 + 𝜌𝑔 23

onde 𝜌 corresponde a massa específica do fluido que atravessa um volume de

controle a uma velocidade �⃗� , 𝜋 corresponde ao tensor tensão e 𝑔 à aceleração da

gravidade.

O problema em estudo consiste em um escoamento incompressível,

axisimétrico (2-D) e sem efeito da gravidade. Desta forma, as equações 22 e 23

podem ser simplificadas respectivamente para as equações 24 e 25, apresentadas

em notação indicial.

0 =

𝜕𝜌𝑢𝑖

𝜕𝑥𝑖 24

𝜌

𝜕𝑢𝑖

𝜕𝑡= −𝜌𝑢𝑖

𝜕𝑢𝑖

𝜕𝑥𝑖+

𝜕𝜋𝑖𝑖

𝜕𝑥𝑖+

𝜕𝜋𝑖𝑗

𝜕𝑥𝑗 25

A análise do escoamento é obtida com a solução das equações 24 e 25 e

com a equação constitutiva para a tensão.

O tensor de tensões 𝜋𝑖𝑖 da equação 25 depende da tensão viscosa,

representado pelo tensor 𝜏, e da pressão p aplicada ao material conforme equação

26.

𝜋 = [

−𝑝 00 −𝑝

] + 𝜏 26

37

Já o tensor de tensão viscosa 𝜏 pode ser especificado pela relação entre a

viscosidade aparente 𝜂 e a taxa de deformação conforme equação 27 para fluido

não newtoniano puramente viscoso [29].

𝜏 = 𝜂 [0

𝜕𝑢𝑥

𝜕𝑟+

𝜕𝑢𝑟

𝜕𝑥𝜕𝑢𝑥

𝜕𝑟+

𝜕𝑢𝑟

𝜕𝑥0

] 27

A viscosidade aparente 𝜂 é calculada conforme equação 20 (item 2.1.3.3 )

para o modelo de tixotropia proposto por de Souza Mendes e cols. [9]. A partir

das simplificações utilizadas e da definição do tensor de tensão 𝜋 (equação 25), da

tensão viscosa 𝜏 (equação 27) e da equação constitutiva para a viscosidade

(equação 20) é possível obter a relação abaixo (equação 28) para a conservação de

quantidade de movimento em função do parâmetro tixotrópico 𝜆.

𝜌

𝜕𝑢𝑖

𝜕𝑡+ 𝜌𝑢𝑖

𝜕𝑢𝑖

𝜕𝑥𝑖= −

𝜕𝑝

𝜕𝑥𝑖+

𝜕

𝜕𝑥𝑗[(

𝜂0

𝜂∞)𝜆

𝜂∞

𝜕𝑢𝑖

𝜕𝑥𝑗] 28

O comportamento de 𝜆 segue a equação constitutiva apresentada pela

equação 21 que depende do valor de 𝜆 em equilíbrio (𝜆𝑒𝑞), além dos parâmetros

característicos do fluido (coeficientes 𝑡𝑒𝑞, 𝑎 e 𝑏). A solução de 𝜆𝑒𝑞, necessária

para resolver a equação constitutiva da tixotropia, pode ser obtida a partir da

relação representada pela equação 29, onde 𝜂𝑒𝑞 é a viscosidade do fluido em

equilíbrio.

𝜆𝑒𝑞 = (ln 𝜂𝑒𝑞 − ln 𝜂∞)/(ln 𝜂0 − ln 𝜂∞) 29

A viscosidade em equilíbrio é obtida experimentalmente da curva de

escoamento. Esse trabalho assume que essa viscosidade pode ser representada

pelo modelo regularizado proposto por de Souza Mendes e cols. [10],

substituindo 𝜂 por 𝜂𝑒𝑞 conforme equação 30.

𝜂𝑒𝑞 = (𝑚�̇�𝑛−1 +

𝜏0

�̇�) (1 − 𝑒

−η0�̇�𝜏0

⁄ ) + 𝜂∞ 30

A Figura 16 resume as equações que devem ser resolvidas e a relação entre

elas, onde fica clara a necessidade em utilizar um método iterativo. Cada uma

destas equações devem ser discretizadas para resolução a partir da ferramenta de

CFD Fluent©.

38

Figura 16 Resumo com as equações que devem ser resolvidas para o escoamento do fluido

tixotrópico

3.2 Discretização da Equação de Transporte Genérica

A discretização das equações diferenciais de transporte que regem o

escoamento apresentadas no item 3.1 , Figura 16, é necessária para a simulação

numérica. Nas ferramentas de CFD, e mais especificamente no Fluent©, o método

utilizado para discretização é o método de volumes finitos. O método de volume

finito consiste resumidamente em: (i) integrar as equações diferenciais em um

volume de controle; (ii) discretizar estas equações; (iii) resolver as equações

linearizadas por um método iterativo, até o resíduo convergir para um valor

mínimo, em cada passo de tempo. [30]

A equação 31 é uma forma genérica para as equações de transporte que

regem o escoamento, e será utilizada como base para apresentar alguns dos

métodos de discretização que podem ser aplicados.

31

onde 𝜌 é a massa específica, 𝜙 é a variável dependente, Γ é o coeficiente difusivo

e 𝑆𝜙 é o termo fonte.

Na equação 31, o termo convectivo expressa o transporte da variável 𝜙

devido ao escoamento e o termo difusivo expressa o transporte por movimento

interno do material, ou seja, por um processo de difusão. A equação de

Eq. 24

Eq. 28

Eq. 20

Eq. 29

Eq. 30

Eq. 24

Eq. 28

Eq. 20

Eq. 29

Eq. 30

𝜕𝜌𝜙

𝜕𝑡+

𝜕𝜌𝜙𝑢𝑖

𝜕𝑥𝑖−

𝜕

𝜕𝑥𝑖(Γ

𝜕𝜙

𝜕𝑥𝑖) = 𝑆𝜙

Termo transiente

Termo advectivo

Termo difusivo

Termo fonte

39

conservação de massa (equação 24) é obtida com ϕ = 1, Sϕ e Γ = 0. Já a

conservação de quantidade de movimento (equação 28) é representada com ϕ =

𝑢𝑖, 𝑆𝜙 = −𝑑𝑝/𝑑𝑥𝑖 e Γ = 𝜂𝑣(𝜆) = (𝜂0

𝜂∞⁄ )

𝜆𝜂∞ .

A integral da equação 31 em um volume de controle genérico (V), aplicando

o teorema da divergência de Gauss, resulta na equação

∫

𝜕𝜌𝜙

𝜕𝑡𝑉

𝑑𝑉 + ∫ 𝜌𝜙𝑢𝑖𝑛𝑖 − Γ𝜕𝜙

𝜕𝑥𝑖𝑛𝑖

𝑠

𝑑𝑆 = ∫ 𝑆𝜙𝑉

𝑑𝑉 32

onde 𝑛𝑖 corresponde a componente normal em cada uma das faces.

A discretização de cada um dos termos da equação 32 será apresentada

resumidamente a seguir.

Diversos métodos podem ser aplicados para discretizar o termo transiente da

equação 32. Um método de primeira ordem pode consistir em aproximar o termo

transiente por uma função 𝐹(𝜙) resultando na igualdade

∫

𝜕𝜌𝜙

𝜕𝑡𝑉

𝑑𝑉 ≈ 𝐹(𝜙) = 𝑉𝜌𝜙𝑛 − 𝜙𝑛−1

∆𝑡 33

onde ∆𝑡 corresponde ao passo de tempo utilizado, 𝜙𝑛 ao valor atual do

termo dependente e 𝜙𝑛−1 ao valor obtido para o tempo anterior.

Já a integral do termo fonte da equação 32 pode ser resolvida a partir do

valor médio em um dado volume, representado por S𝜙̅̅ ̅, conforme equação 34.

∫ S𝜙𝑉

𝑑𝑉 = S𝜙̅̅ ̅𝑉 34

O valor de S𝜙̅̅ ̅ pode ser obtido com a equação de primeiro grau

S𝜙̅̅ ̅ = S𝑐 + S𝑝 𝜙 35

onde S𝑐 e S𝑝 são constantes [30].

O volume de controle genérico em destaque na Figura 17 é composto pelas

faces l, s, n e o, de comprimentos Al, As, An e Ao, respectivamente. O termo

advectivo pode ser aproximado pela equação 36 e o termo difusivo pela equação

37, com base neste volume de controle, a partir do balanço de fluxo que entra e sai

em cada uma das faces.

40

Figura 17 Cálculo da variável nas faces

∫ 𝜌𝜙𝑢𝑖𝑛𝑖𝑠

𝑑𝑆 = 𝜌(𝜙𝑙𝑢𝑙𝐴𝑙 − 𝜙𝑜𝑢𝑜𝐴𝑜 + 𝜙𝑛𝑢𝑛𝐴𝑛 − 𝜙𝑠𝑢𝑠𝐴𝑠) 36

∫ Γ𝜕𝜙

𝜕𝑥𝑖𝑛𝑖

𝑠

𝑑𝑆 = Γ𝑑𝜙

𝑑𝑥|𝑙𝐴𝑙 − Γ

𝑑𝜙

𝑑𝑥|𝑜𝐴𝑜 + Γ

𝑑𝜙

𝑑𝑥|𝑛𝐴𝑛 − Γ

𝑑𝜙

𝑑𝑥|𝑠𝐴𝑠 37

Os valores nas faces necessários para resolver as equações anteriores são

obtidos por interpolação, já que os valores das variáveis são calculados e

armazenados apenas nos centros de cada célula. Diferentes métodos podem ser

aplicados para obter estes valores.

3.2.1 Funções de Interpolação para o Termo Advectivo

No método de primeira ordem, o valor de 𝜙 na face é igual ao valor no

centro da célula a montante desta face. Desta forma, para o exemplo da Figura 18,

o valor de 𝜙𝑙, referente a face entre as células A e B será igual a 𝜙𝐴. Este método

adiciona difusão numérica que pode impactar na convergência requerida e nos

resultados obtidos.

Figura 18 Interpolação no método de primeira ordem para obter o valor da variável 𝝓𝒍,

referente a face entre as células A e B

As

Al Ao

An

Variável φ

Valor de Φ

Sentido do Fluxo

ΦA

ΦB

Cé lula A

Φl

Cé lula B

41

O método de diferenças centrais, mais preciso que o de primeira ordem,

interpola por uma reta os valores da célula a montante e a jusante. Este método

pode resultar em divergência ou oscilações para números de Peclet maiores que

dois [31]. Outro método utilizado é o método power-law que consiste em igualar o

termo convectivo ao difusivo e integrar essa igualdade ao longo do comprimento

Al da face, para escoamento 2-D, obtendo a relação apresentada na equação 38.

𝜙(𝑥) − 𝜙𝑥=0

𝜙𝑥=𝐿 − 𝜙𝑥=0=

𝑒𝑃𝑒𝑥𝐿 − 1

𝑒𝑃𝑒 − 1

38

onde Pe é o número de Peclet, dado por 𝜌𝑢𝐿 Γ𝜙⁄ .

A Figura 19 ilustra o valor de 𝜙 para a face entre duas células para os três

métodos de interpolação descritos acima: primeira ordem, diferenças centrais e

power-law.

Figura 19 Interpolação no método de primeira ordem, diferenças centrais e power-law para

obter o valor da variável 𝝓𝒍, referente a face l entre as células A e B

No método de segunda ordem, o valor da variável na face é obtido a partir

da interpolação dos valores centrais das duas células a montante, em relação ao

sentido de escoamento. Já no método QUICK, a interpolação é realizada com base

nas duas células a montante e na célula a jusante. A Figura 20 ilustra estes dois

métodos.

Valor de Φ

Sentido do Fluxo

ΦA

ΦB

Cé lula A

1

2

3

1: power-law 2: diferenças centrais 3: 1º ordem

Cé lula B

42

Figura 20 Interpolação no método de segunda ordem e QUICK para obter o valor da

variável 𝝓𝒍, referente a face entre as células A e B

3.2.2 Acomplamento da Pressão e Velocidade

Observa-se que o gradiente de pressão precisa ser calculado para resolver a

equação de quantidade de movimento (equação 28). Entretanto, no escoamento

incompressível não é possível calcular diretamente a pressão a partir das equações

fundamentais apresentadas [30]. Desta forma, é necessário resolver indiretamente

a pressão com base nas equações de quantidade de movimento e de continuidade.

Os algoritmos SIMPLE, SIMPLEC ou PISO podem, por exemplo, serem

utilizados para resolver essa questão.

Com o objetivo de controlar a taxa de variação de uma variável entre um

passo de tempo e outro, para aumentar ou evitar divergências na resolução

principalmente de problemas não lineares [8], utiliza-se o conceito de relaxamento

explícito. O novo valor de uma variável 𝜙 depende, desta forma, do valor anterior

desta variável e da variação estimada desta variável, calculada com base em um

dos métodos anteriores, normalizada por um fator ∝ conforme equação 39.

𝜙 = 𝜙𝑎𝑛𝑡𝑒𝑟𝑖𝑜𝑟+∝ ∆𝜙

39

3.2.3 Discretização da Equação de Evolução da Tixotropia

A equação de evolução utilizada neste trabalho segue a forma proposta por

de Souza Mendes e cols. [32], apresentada pela equação 21. Esta equação para ser

utilizada pela ferramenta de CFD deve seguir a forma da equação genérica de

Valor de Φ

Sentido do Fluxo

Cé lula A

2 1: QUICK 2: 2º ordem

Cé lula C

ΦA

ΦB

ΦC

1

Cé lula B

43

transporte dada pela equação 31. Desta forma, o programa poderá discretizar e

resolver a expressão de evolução junto com a solução das demais equações.

Multiplicando ambos os lados da equação 21 pela massa específica

(constante neste caso) e considerando que a derivada no lado esquerdo da equação

equivale a derivada material [25], obtêm-se a equação 40.

40

O termo transiente e advectivo da equação 40 segue a forma da equação

genérica de transporte (equação 31), com termo difusivo nulo, e podem ser

facilmente discretizados e calculados pela ferramenta de CFD. Entretanto, o termo

fonte não segue a relação linear necessária dada pela equação 35. Assim, o termo

fonte desta equação precisa ser linearizado.

A expansão por série de Taylor com truncamento no termo de primeira

ordem parece ser o método mais recomendável para discretização conforme

proposto por Patankar [30]. Este método consiste em estimar o valor do termo de

fonte atual com base nos termos de ordem zero e de primeira ordem da expansão,

ambos referentes aos valores calculados na iteração anterior, conforme equação 41

e ilustrado na Figura 21.

𝑆 = 𝑆∗ +

𝑑𝑆

𝑑𝜆|∗

(𝜆 − 𝜆∗)

41

onde 𝑆∗, 𝑑𝑆 𝑑𝜆⁄ |∗ e 𝜆∗ equivalem respectivamente ao termo fonte, derivada

do termo fonte e parâmetro tixotrópico, todos avaliados na iteração anterior.

Figura 21 Proposta de linearização para o termo fonte. Adaptado de [30]

𝜌𝜕𝜆

𝜕𝑡+ 𝜌

𝜕𝜆𝑢𝑖

𝜕𝑥𝑖 =

𝜌

𝑡𝑒𝑞[(1 − 𝜆)𝑎 − (1 − 𝜆𝑒𝑞)

𝑎(

𝜆

𝜆𝑒𝑞)

𝑏

]

Termo transiente

Termo advectivo

Termo fonte

44

Aplicando a equação 41 ao termo fonte da equação 40, obtêm-se a igualdade

representada pela equação abaixo.

𝑆 =

𝜌

𝑡𝑒𝑞[𝑎(𝜆 − 1)

(𝜆∗ − 1)+ (1 − 𝑎)] [−(𝜆∗ − 1)]𝑎

+𝜌

𝑡𝑒𝑞(−(𝜆∗

𝑒𝑞 − 1))𝑎[(𝑏 − 1) −

𝑏𝜆

𝜆∗ ] (𝜆∗

𝜆∗𝑒𝑞

)

𝑏

42

Separando os termos da equação 42 é possível obter os valores de 𝑆𝑐 e

𝑆𝑝 da equação 35. As equações 43 e 44 representam, respectivamente, os termos

𝑆𝑝 e 𝑆𝑐 obtidos.

S𝑝 =𝜌

𝑡𝑒𝑞

𝑎𝜆∗[−(𝜆∗ − 1)]𝑎 − 𝑏(𝜆∗ − 1)[−(𝜆∗𝑒𝑞 − 1)]

𝑎(

𝜆∗

𝜆∗𝑒𝑞)𝑏

𝜆∗(𝜆∗ − 1)

43

S𝑐 = −

𝜌

𝑡𝑒𝑞(𝜆∗ − 1)

{[(𝑎 − 1)𝜆∗ + 1][−(𝜆∗ − 1)]𝑎

− (𝑏 − 1)(𝜆∗ − 1)[−(𝜆∗𝑒𝑞 − 1)]

𝑎(

𝜆∗

𝜆∗𝑒𝑞

)

𝑏

}

44

Com os termos S𝑝 e S𝑐 calculados têm-se a equação evolutiva de tixotropia

escolhida de acordo com a forma e premissas da equação genérica de transporte

(equação 31), que pode ser resolvida pela ferramenta de CFD.

45

4 Metodologia

4.1 Algoritmo para Resolução das Equações de Transporte no Fluent©

A ferramenta de CFD Fluent© é capaz de resolver as equações de transporte

discretizadas utilizando dois métodos numéricos, a saber: por pressão ou por

densidade. No método por densidade, a equação da continuidade é resolvida junto

com a equação de quantidade de movimento e o campo de pressão é determinado

a partir de uma equação de estado. Já no método por pressão, a pressão é obtida a

partir das equações de continuidade e quantidade de movimento [8].

O método por pressão, que será utilizado para resolver o problema aqui

proposto, pode ser resolvido utilizando o algoritmo acoplado ou segregado. No

algoritmo acoplado as equações de conservação de quantidade de movimento e

continuidade são resolvidas simultaneamente. No método segregado, primeiro as

equações de quantidade de movimento são resolvidas e, em seguida, a de

continuidade e demais, tais como a de conservação de energia, turbulência ou

aquelas definidas pelo usuário. A Figura 22 apresenta a sequência de operações

para o algoritmo segregado e acoplado. Em uma simulação transiente, esta

sequência é aplicada para cada passo de tempo.

Figura 22 Algoritmo para resolução das equações de transporte no Fluent©. Adaptado de [8]

Verifica Convergência

Atualiza as Propriedades

Resolve Equações de Movimento em x, y e z

Resolve Equação da Continuidade

Atualiza Fluxo de Massa, Velocidades e Pressão

Resolve Equações de Energia, Turbulência, etc.

Convergiu

Algoritmo Segregado

Não

Convergiu

Atualiza as Propriedades

Resolve Equações de Movimento e Continuidade

Atualiza Fluxo de Massa

Resolve Equações de Energia, Turbulência, etc.

Verifica Convergência

Convergiu

Algoritmo Acoplado

Não

Convergiu

46

Apesar de o algoritmo acoplado apresentar uma taxa de convergência mais

elevada [8], o algoritmo segregado foi utilizado já que nos primeiros testes

realizados apresentou soluções mais rápidas.

O simulador Fluent© não é capaz de simular nativamente o escoamento

para fluidos tixotrópicos. Entretanto, é possível via funções escritas pelo usuário

inserir a resolução de uma equação de transporte qualquer, desde que siga a forma

apresentada no item 3.2 . A partir desta funcionalidade, a equação de transporte

para a tixotropia pode ser programada e acoplada à equação de conservação de

quantidade de movimento com a equação de viscosidade apropriada. Algumas

funções já existentes na biblioteca do Fluent© permitem este acoplamento. A

Figura 23 apresenta as funções utilizadas para permitir a simulação numérica do

escoamento de um fluido com comportamento tixotrópico.

Estas funções, escritas em linguagem C, podem ser compiladas pelo próprio

Fluent© antes de iniciar a simulação ou podem ser interpretadas durante a

simulação. Compilar as funções reduz o tempo de execução da simulação e

permitem utilizar funções mais complexas e necessárias para simular a tixotropia.

Figura 23 Funções adicionais requeridas para implementar a tixotropia no Fluent© para o

algoritmo segregado

O termo transiente, advectivo e fonte da equação de transporte para a

tixotropia (equação 40) são descritas em linguagem C utilizando respectivamente

Verifica Convergência

Atualiza as Propriedades

Resolve Equações de Movimento em x, y e z

Resolve Equação da Continuidade

Atualiza Fluxo de Massa, Velocidades e Pressão

Resolve Equações de Energia, Turbulência, etc.

Convergiu

Algoritmo Segregado

Não

Convergiu

Equação de Evolução da Tixotropia

Macro SOURCE

Macro UNSTEADY

Macro FLUX

Macro PROPERTY Macro ADJUST

47

as funções Unsteady, Flux e Source. Os algoritmos utilizados são apresentados

nas tabelas a seguir.

Tabela 1 Algoritmo para o termo fonte

Tabela 2 Algoritmo para o termo de fluxo

Tabela 3 Algoritmo para o termo transiente

A relação entre a viscosidade aparente e o termo tixotrópico, equação 20

(item 2.1.3.3), é descrita na função Property conforme lógica da Tabela 4.

Tabela 4 Algoritmo para a viscosidade

O termo fonte que foi linearizado conforme metodologia apresentada no

item 3.2.3 depende do valor da iteração anterior. A função Adjust, que é executada

ao final dos cálculos de uma iteração, foi necessária para guardar o valor atual de

algumas variáveis para a próxima iteração. O algoritmo é apresentado na Tabela

5.

Tabela 5 Algoritmo para manipular valores anteriores

DEFINE_SOURCE 1. Obter taxa de cisalhamento na célula com a função

C_STRAIN_RATE_MAG; 2. Obter valores da iteração anterior em C_UDMI 3. Calcular viscosidade de equilíbrio; 4. Calcular estrutura de equilíbrio; 5. Calcular o termo fonte; 6. Salvar derivada do termo fonte; 7. Retornar termo fonte calculado.

DEFINE_UDS_FLUX 1. Utilizar macro F_FLUX (resolve 𝜌�⃗� ) para retornar fluxo.

DEFINE_UDS_UNSTEADY 1. Obter passo de tempo atual com função RP_Get_Real; 2. Obter volume da célula atual com a função C_VOLUME; 3. Obter a massa na célula atual com a função C_R; 4. Calcular os coeficientes do termo fonte e central da

equação diferencial discretizada.

DEFINE_PROPERTY 1. Obter o valor da estrutura da célula salvo em C_UDSI; 2. Calcular e retornar o valor da viscosidade.

DEFINE_ADJUST 1. Obter ponteiros para thread (via Thread *) e célula (via

cell_t); 2. Repetir para cada thread e cada célula:

i. Salvar valor das variáveis atuais na memória com C_UDMI;

ii. Atualizar com valores das variáveis anteriores.

48

4.2 Principais Configurações no Fluent©

A simulação numérica utilizou a versão 16.0.0 do Fluent© com a opção de

precisão dupla habilitada. O método de pressão foi utilizado para resolver as

equações de transporte (pressure-based solver) para o modelo transiente

axisimétrico com geometria de duas dimensões (2-D) e as seguintes

configurações:

Condição inicial: fluido parado no duto e totalmente estruturado (λ=1).

Condições de entrada:

- Degrau de pressão aplicado no instante inicial, pressão mantida

constante em toda a simulação;

- Gradiente do escalar λ igual a zero na direção normal a entrada

conforme utilizado por Wachs e cols. [23].

Condições de saída

- Pressão constante e igual a zero;

- Gradiente do escalar λ igual a zero na direção normal a saída.

Condições na parede

- Sem escorregamento na parede (slip igual a zero);

Passo de tempo: igual a 0,0001s.

O algoritmo PISO para acoplamento da pressão e velocidade foi escolhido

por ser o mais recomendável em simulações transientes e oferecer estabilidade

mesmo para relaxamentos elevados [8]. Os fatores de relaxamento apresentados

na Tabela 6 foram utilizados e obtidos após algumas simulações preliminares

indicarem menor tempo de convergência e estabilidade para estes fatores.

Tabela 6 Fatores de relaxação utilizados

Variável Fator de

relaxamento

Pressão 0.4

Densidade 1

Forças de corpo 1

Momento 0.4

Escalar 1

Os métodos de solução utilizados para discretização, conforme descrito

brevemente no item 3.2, “Discretização da Equação de Transporte Genérica”,

foram:

49

Cálculo de pressão: 2ª ordem;

Cálculo do momento: Power-law;

Cálculo do escalar (parâmetro tixotrópico): 2ª ordem;

Cálculo dos gradientes: mínimos quadrados.

Os critérios de convergência utilizados para cada um dos resíduos das

equações de conservação de massa, momento e tixotropia são apresentados na

Tabela 7.

Tabela 7 Resíduos máximos (em valor absoluto)

Equação Resíduo

Máximo

Continuidade 1.10-3

Velocidade Axial 1.10-5

Velocidade Radial 1.10-5

Escalar 1.10-5

4.3 Teste de Independência de Malha

O impacto da malha utilizada no comportamento da simulação e nos

resultados gerados foram avaliados em regime permanente para alguns cenários

utilizando a equação 10 da tixotropia (item 2.1.1.5) com a e b iguais a 1. A

geometria de teste representa um duto de comprimento 2 m e raio de 0.05 m,

resultando em uma razão comprimento por largura de 40. A Tabela 8 apresenta a

quantidade de células em cada uma das malhas utilizadas.

Tabela 8 Diferentes malhas utilizadas

Nome da Malha

Número de Células

Comprimento Altura Total

1 100 5 500

2 200 10 2000

3 400 20 8000

4 640 32 20480

5 720 36 25920

6 1080 54 58320

7 1440 72 103680

8 2160 180 388800

Especificamente para esta análise, foi utilizada a convergência mínima de

1.10-7 para as equações de momento, continuidade e parâmetro construtivo. Os

valores de velocidade total média a 10 cm da saída (L=1.90 m) são apresentados

para cada uma das malhas utilizadas na Figura 24. A velocidade média calculada

50

para a malha mais refinada (2160x180) é aproximadamente igual aquela com

refinamento imediatamente abaixo (1440x72).

Figura 24 Velocidade Axial média a 10cm da saída em função do número de células

utilizadas na rede

A Figura 25 apresenta o erro relativo à malha mais fina para a velocidade

média (linha vermelha) e para a velocidade em cada célula ao longo do raio (linha

azul). A malha 5, de 720 por 36 células, foi escolhida. Além de possuir menos que

a metade do número de células da malha 6, resultou em erro relativo de

aproximadamente 2% em relação à malha mais fina (contra 1% para a malha 6).

Figura 25 Erros relativos à malha mais refinada (2160x180).

0,93

0,94

0,95

0,96

0,97

0,98

0,99

1,00

1,01

100 1000 10000 100000 1000000

Ve

l. M

éd

ia \

Ve

l. M

éd

ia M

áxim

a

Nº de Células

Velocidade Média pelo Número de Células da Malha

0%

1%

2%

3%

4%

5%

6%

7%

8%

9%

10%

100 1000 10000 100000 1000000

Erro

(%

)

Nº de Células

Erro Relativo para as Velocidades e Velocidades Médias

Erro Médio das Velocidades Erro da Velocidade Média

Malha escolhida

1

2

3

4

5

6

7 8

51

Com o objetivo de validar a malha escolhida foi utilizada a solução analítica

para escoamento de um fluido que segue o modelo de Bingham (item 2.1.1.3 )

para a viscosidade.

O perfil da velocidade 𝑢 para escoamento laminar desenvolvido e em

regime permanente para escoamento de um fluido de Bingham é dado pela

equação

𝑢(𝑟) = {𝑢 =

1

𝑘{(𝑅2 − 𝑟2)

Δ𝑝

4L+ (𝑟 − 𝑅)𝜏0} 𝑝𝑎𝑟𝑎 �̇� > �̇�𝑐

Δ𝑝

4Lk(𝑅 − 𝑟0)

2 𝑝𝑎𝑟𝑎 �̇� = �̇�𝑐

45

onde Δp representa o diferencial de pressão, R o raio do duto, L o comprimento do

duto, k a constante de Bingham e 𝑟0 a posição onde a taxa de deformação é igual a

taxa crítica, i.e., a posição em que o módulo da tensão é igual a tensão limite de

escoamento.

O resultado do cálculo analítico utilizando a equação 45 é comparado com o

resultado no Fluent com o modelo de viscosidade de Bingham implementado. A

Figura 26 apresenta os perfis de velocidade para ambos os resultados e o erro

entre eles, inferior a 2% e considerado adequado.

Figura 26 Perfil de velocidade obtida pelo cálculo analítico e no modelo do Fluent

-2,50%

-2,00%

-1,50%

-1,00%

-0,50%

0,00%

0,50%

1,00%

1,50%

0

0,1

0,2

0,3

0,4

0,5

0,6

0,7

0,8

0 0,01 0,02 0,03 0,04 0,05

Erro

(%

)

Ve

loci

dad

e (

m/s

)

Raio (m)

Analítico

Fluent

Erro

52

5 Resultados e Discussão

5.1 Propriedades Reológicas do Fluido Utilizado

O fluido utilizado como referência nas simulações foi uma solução de

Laponita de concentração 0,005 mol/l de sal, conforme discutido por Lima [3]. O

modelo regularizado de viscosidade proposto por de Souza Mendes [10], equação

10 do Item 2.1.1.5, foi utilizado. Os parâmetros reológicos para este modelo (vide

Tabela 9) são os mesmos que os apresentados por Lima [3] para esta solução de

Laponita.

Tabela 9 Parâmetros Reológicos da Solução de Laponita Utilizada como Referência

Parâmetro Valor

𝑚 1,9 Pa

𝑛 0,2717

𝜏0 4,6 Pa

𝜂0 1,0E+6 Pa.s

𝜂∞ 0,001 Pa.s

Os valores utilizados para a equação de evolução da tixotropia adotada, a, b

e 𝑡𝑒𝑞 da equação 21, foram definidos arbitrariamente por falta de dados

experimentais. Para o caso base, os valores de a e b foram 0,5 e o valor para o

tempo de equilíbrio (𝑡𝑒𝑞) foi igual a 0,6. A Figura 27 a seguir compara a taxa de

variação da estrutura (λ) do fluido com o tempo para o caso base (a e b iguais a

0,5) e para o caso linear (a e b iguais a 1,0).

Figura 27 Variação da estrutura λ do fluido para λeq igual a 0,5 em função dos parâmetros a

e b

0,00

0,10

0,20

0,30

0,40

0,50

0,60

0,70

0,80

0,90

1,00

0,00 0,20 0,40 0,60 0,80 1,00

dλ/

dt

(1/s

)

λ

dλ/dt

linear / λeq=0.5

a=b=0.5 / λeq=0.5

53

5.2 Análises Qualitativas para o Reinício

O algoritmo apresentado no item 4.1 será testado para um cenário típico de

reinício de escoamento. O fluido se encontra em repouso e totalmente estrutura

conforme condições apresentadas no item 4.2 e, no instante inicial é aplicado um

degrau de pressão na entrada, e a velocidade média na saída próximo a saída é

medido, conforme ilustra a Figura 30.

Figura 28 Cenário inicial apresentando a pressão aplicada na entrada e a área rachurada do

duto em que a simulação é realizada (sistema axisimétrico)

Dois resultados serão apresentados, com dois valores de pressão aplicada na

entrada. No cenário de baixa pressão, a condição de contorno adotada neste

trabalho, gradiente de λ igual a zero conforme item 4.2 , é comparada com a

condição de contorno com o valor de λ igual a 0,2.

5.2.1 Pressão de Entrada Baixa

A evolução da velocidade média obtida para um degrau de pressão aplicado

no instante inicial é apresentada na figura abaixo. Observa-se que o sistema

demora a responder e a velocidade aparente do fluido aumenta somente após 6s.

Figura 29 Velocidade Média do Fluido

0,00

0,02

0,04

0,06

0,08

0,10

0,12

0,14

0,16

0,00 10,00 20,00 30,00 40,00 50,00 60,00

Ve

loci

dad

e M

éd

ia (

m/s

)

Tempo (s)

Velocidade Média do Fluido

Velocidade (x/L=0.975)

54

A Figura 30 apresenta o mapa da estrutura do fluido (representada pelo

escalar λ) desde o instante inicial, com fluido aproximadamente em repouso, até o

sistema atingir regime permanente, após a aplicação de uma pressão baixa na

entrada (1500 Pa).

Figura 30 Mapa da estrutura do fluido para diferentes tempos ao longo do reinício com

pressão de entrada de 1500 Pa e teq de 0,6 s

No instante inicial, o sistema está em repouso e a estrutura do fluido é

máxima e igual a um (Figura 30A). Na teoria este fluido se comportaria como

sólido com viscosidade infinita. Entretanto, devido a limitações da simulação

considerada, é necessário simplificar este comportamento de sólido

representando-o por um fluido de alta viscosidade (no caso igual a 1.106 Pa.s)

conforme apresentado no item 2.1.1.5.

Imediatamente é aplicada uma pressão na entrada igual a 1500 Pa. Neste

momento, o fluido começa a se mover a baixíssima velocidade (devido à

representação simplificada de sólido via alta viscosidade). Este movimento causa

cisalhamento do fluido com a parede e inicia-se a quebra de sua estrutura próxima

A (t=0.5s)

B (t=10.5s)

C (t=20.5s)

D (t=30.5s)

E (t=50.5s)

F (t=60.5s)

55

à parede (Figura 30B). A estrutura do fluido nos primeiros instantes de tempo é

apresentada na Figura 31 abaixo.

Figura 31 Perfil no raio para a estrutura do fluido em diferentes tempos com pressão de

entrada de 1500 Pa e teq de 0,6 s

A velocidade do escoamento já é significativa após alguns segundos, se

comparada à velocidade final, e o fluido continua a se quebrar próximo à parede

com o aumento gradual da velocidade (Figura 30C e D). Observa-se também a

formação de um vale no parâmetro λ, que se propaga com o tempo, causado

devido às perturbações na taxa de deformação próximo à entrada (escoamento não

desenvolvido). Observa-se na Figura 32 que a taxa de deformação no centro sofre

aumento significativo a partir de 14 s na região próxima à entrada.

Figura 32 Taxa de deformação no centro em diferentes posições ao longo do comprimento

0,20

0,40

0,60

0,80

1,00

0,00 2,00 4,00 6,00 8,00 10,00 12,00

Estr

utu

ra d

o F

luid

o

Tempo (s)

Estrutura do Fluido em x=3.90 (x/L=0.975)

Centro (0/5)

r=0.01m (1/5)

r=0.02m (2/5)

r=0.03m (3/5)

r=0.04m (4/5)

Parede (5/5)

1,0E-08

1,0E-07

1,0E-06

1,0E-05

1,0E-04

1,0E-03

1,0E-02

0,00 10,00 20,00 30,00 40,00 50,00 60,00 70,00

Taxa

de

Cis

alh

ame

nto

(1

/s)

Tempo (s)

Cisalhamento no Centro

x/L = 0.25

x/L = 0.5

x/L = 0.75

x/L = 0.875

x/L = 0.975

56

Com o passar do tempo, todo o fluido que estava inicialmente no duto é

removido (Figura 30E). A estrutura do fluido e a velocidade média entre o tempo

de 50,5 s (Figura 30E) e de 60,5 s (Figura 30F) sofre pouca alteração e tende ao

valor constante conforme Figura 33: o sistema atingiu regime permanente.

Figura 33 Estrutura do fluido no centro em diferentes instantes de tempo

5.2.1.1 Resultados para Valor de λ Constante na Entrada

Neste item será avaliado o resultado para um cenário com condição de

contorno de λ igual a 0,2. O objetivo é comparar este resultado com os valores

obtidos para o caso do grandiente de λ igual a zero.

A velocidade obtida para a simulação de reinício com λ na entrada igual a

0,2 é apresentada na figura abaixo

,

junto com a velocidade obtida para o cenário anterior (baixa pressão com

0,00

0,05

0,10

0,15

0,50

0,60

0,70

0,80

0,90

1,00

0,00 10,00 20,00 30,00 40,00 50,00 60,00

Ve

loci

dad

e M

éd

ia (

m/s

)

Estr

utu

ra M

éd

ia λ

Tempo (s)

Velocidade Média e Estrutura Média do Fluido

λ x/L=0.25 λ x/L=0.50

λ x/L=0.75 λ x/L=0.875

λ x/L=0.975 Velocidade (x/L=0.975)

0,00

0,02

0,04

0,06

0,08

0,10

0,12

0,14

0,16

0,00 10,00 20,00 30,00 40,00 50,00

Ve

loci

dad

e M

éd

a (m

/s)

Tempo (s)

Velocidade Média do Fluido para as Duas Condições de Contorno de λ

Gradiente igual a zero -Velocidade (x/L=0.975)

Valor constante igual a 0,2- Velocidade (x/L=0.975)

57

grandiente de λ igual a zero na entrada, conforme item 5.2.1 ). O valor da pressão

aplicada neste item foi ajustado para que a velocidade final em regime permanente

se aproximasse do cenário anterior. Observa-se que qualitativamente o

comportamento da velocidade final não sofre influência da condição da entrada.

Figura 34 Velocidade Média do Fluido para gradiente de λ igual a zero e λ igual a 0,2 como

condição de contorno na entrada

Verifica-se a partir da Figura 35 que o comportamento da estrutura do fluido

com o tempo segue a mesma tendência daquele observado para o cenário anterior,

com gradiente de λ igual a 0.

0,00

0,02

0,04

0,06

0,08

0,10

0,12

0,14

0,16

0,00 10,00 20,00 30,00 40,00 50,00

Ve

loci

dad

e M

éd

a (m

/s)

Tempo (s)

Velocidade Média do Fluido para as Duas Condições de Contorno de λ

Gradiente igual a zero -Velocidade (x/L=0.975)

Valor constante igual a 0,2- Velocidade (x/L=0.975)

58

Figura 35 Mapa da estrutura do fluido para diferentes tempos ao longo do reinício com

condição de contorno na entrada para λ igual a 0,2 e teq de 0,6 s

Os distúrbios encontrados a partir do tempo de 20,5s para o cenário com

gradiente de λ igual a zero na entrada (Figura 30C) não são observados na

simulação com igual λ a 0,2 Figura 35C). Na região próximo à saída, onde o

sistema atinge o escoamento desenvolvido, o resultado é semelhante para os dois

casos.

Observa-se que o fluido entrante desloca aquele inicialmente em repouso ao

longo do tempo e, enquanto caminha pelo duto, sua estrutura se regenera até

atingir os valores verificados para o sistema desenvolvido (trecho final do duto,

Figura 38E e F). A Figura 36 detalha a condição de λ com o tempo em diferentes

posições e ilustra o efeito do fluido entrante no sistema, à semelhança do cenário

com gradiente de igual a zero na entrada (Figura 33).

A (t=0.5s)

B (t=10.5s)

C (t=20.5s)

D (t=30.5s)

E (t=50.5s)

F (t=60.5s)

59

Figura 36 Estrutura do fluido no centro em diferentes instantes de tempo para λ igual a 0,2

na entrada

5.2.2 Pressão de Entrada Alta

A mesma condição anterior, repartida de um fluido em repouso, é simulada

neste item para uma pressão maior de entrada (4500 Pa). Conforme discutido no

item 5.2.2, as velocidades obtidas para o regime permanente ultrapassam o limite

de escoamento laminar. Entretanto, até 6 segundos depois de aplicada a pressão,

quando a velocidade ainda é pequena, o fluxo é laminar conforme Figura 37.

Figura 37 Número de Reynolds e Comprimento de Entrada calculados próximo a saída ao

longo do tempo para pressão aplicada alta (4500Pa)

0,00

0,05

0,10

0,15

0,20

0,50

0,60

0,70

0,80

0,90

1,00

0,00 10,00 20,00 30,00 40,00 50,00 60,00

Ve

loci

dad

e M

éd

ia (

m/s

)

Estr

utu

ra M

éd

ia λ

Tempo (s)

Velocidade Média e Estrutura Média do Fluido

λ x/L=0.25 λ x/L=0.50

λ x/L=0.75 λ x/L=0.875

λ x/L=0.975 Velocidade (x/L=0.975)

0

5

10

15

20

25

30

0

1000

2000

3000

4000

5000

6000

0,00 2,00 4,00 6,00 8,00 10,00 12,00 14,00 16,00 18,00 Co

mp

rim

en

to d

e E

ntr

ada

(m)

Nú

me

ro d

e R

eyn

old

s

Tempo (s)

Número de Reynolds e Comprimento de Entrada Próximo à Saída (x/L=0.975)

Número de Reynolds Comprimento de Entrada

60

O mapa para a estrutura do fluido ao longo do reinício com pressão de

entrada igual a 4500 Pa é apresentado na Figura 38.

Figura 38 Mapa da estrutura do fluido para diferentes tempos ao longo do reinício com

pressão de entrada de 4500 Pa e teq de 0,6 s

Em comparação ao cenário de baixa pressão (Figura 30, item 5.3.1), os

seguintes pontos são destacados em relação a partida com alta pressão aplicada:

Nos primeiros instantes de tempo (Figura 38A e B), a estrutura do

fluido é reduzida a partir da parede, onde a taxa de deformação é

mais alta. O perfil de velocidade sofre a alteração ilustrada na Figura

39 para os primeiros instantes.

Em um dado momento um vale na estrutura do fluido se forma

(Figura 38C), assim como observado para baixa pressão. Além

disso, devido às velocidades mais elevadas, observa-se oscilação na

taxa de deformação que resulta na oscilação da estrutura do fluido,

conforme trecho final da Figura 38C e Figura 38D.

A (t=0.5s)

B (t=3.0s)

C (t=4.5s)

D (t=6.5s)

E (t=16s)

F (t=26s)

61

Com o passar do tempo o fluido que entra expulsa todo o fluido

original (Figura 38E) e o sistema atinge o regime permanente

(Figura 38F).

Figura 39 Perfil de velocidade em diferentes instantes de tempo para pressão elevada

aplicada

5.3 Avaliação das Configurações Iniciais

5.3.1 Condição de Contorno na Entrada para Lambda

O problema proposto consiste em um trecho de tubulação onde um fluido

tixotrópico totalmente estruturado repousa. A condição de contorno para a

estrutura do fluido na entrada, representada pelo parâmetro λ, pode assumir uma

dentre as três opções a seguir: i) valor para λ constante ou igual a uma função

definida pelo usuário; ii) valor do fluxo de λ igual a uma constante ou função; ou

iii) outro fluido na entrada (não tixotrópico), em uma simulação multifásica. A

terceira opção, simulação multifásica, não será considerada neste momento por

extrapolar os objetivos deste trabalho.

Com o objetivo de avaliar as duas opções para a condição de contorno de λ,

as seguintes análises foram realizadas para um degrau de pressão aplicado de

1500Pa na entrada:

i. Valor da estrutura na entrada constante e igual a 0,2 e igual a 0,4;

ii. Valor do fluxo igual a zero, em semelhança à condição para a

velocidade em um escoamento desenvolvido.

0,00

0,20

0,40

0,60

0,80

1,00

1,20

0,00 0,01 0,02 0,03 0,04 0,05

Ve

loci

dad

e /

Ve

loci

dad

e M

áxim

a

Raio (m)

Velocidade normalizada em diferentes instantes medida na saída (x/L=0.975)

t=1s

t=2s

t=3s

t=4s

t=5s

t=6s

62

A partir da velocidade média do fluido, obtida próxima da saída para evitar

perturbações de entrada, verifica-se que a solução em regime permanente depende

da condição de contorno para a estrutura, conforme Figura 40. Os valores para λ

na entrada de 0,2 e 0,4 foram escolhidos para resultar em velocidades médias

próximas à condição de fluxo igual a zero.

Figura 40 Velocidade média do fluido a 3,5 m da entrada para pressão de entrada de 1500 Pa

e teq de 0,6 s para diferentes condições de entrada para λ

Independente das três condições de entrada para λ consideradas, a solução

para o início do escoamento aparenta ser aproximadamente a mesma. Isso ocorre

porque o fluido novo (fluido com mesma característica que aquele inicialmente da

tubulação, mas com diferente história) que entra na tubulação, afetado pela

condição de λ, ainda não deslocou significativamente o fluido original e, desta

forma, não afetou ainda o escoamento. A solução sofre pouca influência da

condição de contorno de λ até o tempo 18 s, quando o fluido novo ainda não

ocupou 5% do total do duto conforme Figura 41.

Quando o novo fluido já ocupa a maior parte do duto, observa-se que um

menor valor de λ na entrada resulta em maior velocidade final porque o fluido

novo entra com menor viscosidade. O escoamento atinge regime permanente no

tempo final aproximado de 48 segundos para as três simulações.

0,00

0,02

0,04

0,06

0,08

0,10

0,12

0,14

0,16

0,00 10,00 20,00 30,00 40,00 50,00

Ve

loci

dad

e M

éd

ia (

m/s

)

Tempo (s)

Velocidade Média por Condição de Entrada de λ

∇λ = 0

λ = 0.2

λ = 0.4

63

Figura 41 Posição média do fluido a partir da integral da velocidade média a 3,5m da

entrada para diferentes condições de entrada para λ

A velocidade total no centro ao longo do comprimento é apresentada na

Figura 42 após atingir o regime permanente. Quando o fluxo de λ é imposto a

velocidade atinge o regime desenvolvido rapidamente, velocidade no centro não

varia ao longo do comprimento, o que indica ser mais adequada esta condição de

λ do que as outras opções analisadas.

Figura 42 Velocidade do fluido ao longo do comprimento para a pressão de entrada de 1500

Pa e teq de 0,6 s para diferentes condições de entrada para λ no tempo final

O desenvolvimento da estrutura do fluido ao longo do comprimento é mais

suave e tende para o estado permanente para a condição na entrada de fluxo de λ

igual a zero. Quando o valor de λ é definido para a entrada, o parâmetro

tixotrópico não atinge o estado permanente ou sofre uma variação muito brusca

0,00

0,50

1,00

1,50

2,00

2,50

3,00

3,50

4,00

0,00 10,00 20,00 30,00 40,00 50,00

Dis

tân

cia

Mé

dia

Pe

rco

rrid

a (m

)

Tempo (s)

Distância Média Percorrida

∇λ = 0

λ = 0.2

λ = 0.4

0,10

0,12

0,14

0,16

0,18

0,20

0,22

0,24

0,26

0,28

0,00 0,50 1,00 1,50 2,00 2,50 3,00 3,50 4,00

Ve

loci

dad

e (

m/s

)

Comprimento (m)

Velocidade no Centro (raio = 0) por Comprimento

∇λ =0

λ =0.2

λ =0.4

64

próxima à saída conforme ilustra a Figura 43. As condições de contorno para λ na

saída são iguais para os três casos analisados (fluxo igual a zero).

Figura 43 Estrutura do fluido para diferentes condições de entrada de λ no tempo final para

a pressão aplicada de 1500 Pa e teq de 0,6 s

Devido aos resultados anteriores, a condição de contorno escolhida para a

tixotropia na entrada é com fluxo de λ igual a zero. Todos os resultados

apresentados nos itens a seguir consideram esta condição. Conforme discutido,

essa condição não afeta a análise de repartida (primeiros instantes de tempo), mas

altera a evolução do fluido até atingir o regime permanente.

5.3.2 Análise do Comprimento da Malha

As primeiras simulações realizadas consideravam um duto de comprimento

de 2 m e razão comprimento-diâmetro, L/D, igual a 20. Devido à dificuldades em

obter convergência em algumas análises preliminares, a razão L/D foi dobrada.

Entretanto, foi identificado que outras questões no código e na configuração do

simulador eram causadoras deste problema após este aumento de comprimento.

De qualquer forma, o aumento da razão L/D é válido, pois pode permitir que as

perturbações de entrada, que afetam tanto a velocidade como a estrutura do fluido,

sejam suavizadas ao longo do comprimento.

A velocidade média em diferentes posições foi comparada à velocidade

média próxima à saída quando o sistema atinge o regime permanente. As pressões

aplicadas dependem da razão comprimento-diâmetro (L/D) utilizada e foram

0,10

0,20

0,30

0,40

0,50

0,60

0,70

0,80

0,00 0,50 1,00 1,50 2,00 2,50 3,00 3,50 4,00

Estr

utu

ra M

éd

ia d

o F

luid

o

Comprimento (m)

Estrutura Média do Fluido por Comprimento

∇λ =0

λ =0.2

λ =0.4

65

calculadas de tal maneira que resultem em velocidades aproximadamente iguais

como pode ser observado na Figura 44.

Figura 44 Velocidade média no tempo para as duas razões comprimento-diâmetro

(velocidade medida na razão posição comprimento (x/L) igual a 0.975)

O escoamento atinge regime desenvolvido para estes cenários em ambas as

razões L/D analisadas. O resultado é apresentado na Figura 45, onde a diferença

entre a velocidade no centro em diferentes posições e a velocidade no centro

próximo à saída é ilustrada. A variação da velocidade com a posição para a razão

L/D igual a 20 é significativamente maior que com razão L/D igual a 40. Por essa

razão, o comprimento de 4 m (L/D igual a 40) foi adotado.

Figura 45 Velocidade Média em relação a velocidade média final para diferentes razões

posição-comprimento para L/D 40 e L/D 20

0,00

0,02

0,04

0,06

0,08

0,10

0,12

0,14

0,16

0,00 5,00 10,00 15,00 20,00 25,00 30,00 35,00 40,00

Ve

loci

dad

e M

éd

ia (

m/s

)

Tempo (s)

Velocidade Média Próxima a Saída (x/L=0.975)

L/D=40 L/D=20

-0,2%

0,0%

0,2%

0,4%

0,6%

0,8%

1,0%

1,2%

x/L = 0.25 x/L = 0.5 x/L = 0.75 x/L = 0.875 x/L = 0.975

% V

elo

cid

ade

no

po

nto

de

saí

da

(x/L

=0,9

75

)

Razão Posição-Comprimento x/L

Velocidade no Centro em Diferentes Posições para o Tempo Final

L/D=20

L/D=40

66

A taxa de deformação no centro esperada para escoamento desenvolvido é

zero. A taxa de deformação no centro é inferior a 10-6s-1 para L/D=40 (vide Figura

46), que pode ser considerada desprezível, em todo o tempo de simulação. Com

L/D=20 a taxa de deformação atinge valores superiores a 10-5s-1, que também é

baixa, mas uma ordem de grandeza superior aquela com L/D=40.

Figura 46 Taxa de deformação no centro e velocidade média para pressão baixa aplicada e

razão L/D igual a 40 e 20

Algumas simulações foram realizadas com diferentes pressões aplicadas

para comprimento de duto igual a 4m (razão L/D igual a 40). Com o objetivo de

verificar se nestas simulações o regime é laminar e se o comprimento do duto é

suficiente para atingir escoamento desenvolvido, foi utilizado simplificadamente o

número de Reynolds para fluido newtoniano, utilizando a viscosidade na parede

(mínima) e a velocidade média no regime permanente. A partir do número de

Reynolds calculado, observa-se que o regime ultrapassa a faixa laminar (Reynolds

igual a 2300) para a maioria dos casos. A Figura 47 apresenta o número de

Reynolds e o comprimento de entrada para estas simulações.

0,00

0,02

0,04

0,06

0,08

0,10

0,12

0,14

0,16

1,0E-08

1,0E-07

1,0E-06

1,0E-05

1,0E-04

0,00 10,00 20,00 30,00 40,00 50,00 60,00

Ve

loci

dad

e M

éd

ia (

m/s

)

Taxa

de

Cis

alh

ame

nto

(1

/s)

Tempo (s)

Pressão Baixa - Cisalhamento no Centro e Velocidade Média Próxima a Saída (x/L=0.975)

Taxa (L/D=40) Taxa (L/D=20)

Velocidade (L/D=40) Velocidade (L/D=20)

67

Figura 47 Número de Reynolds e Comprimento de Entrada, para diferentes pressões

aplicadas, calculados próximo a saída

Somente três cenários de pressão aplicada iguais a 1500, 2000 e 2500 Pa

resultam em escoamento laminar considerando o número de Reynolds calculado.

Desta forma, apenas a análise em regime permanente destes cenários deve ser

válida. A Figura 47 mostra também que o comprimento de entrada a partir da

pressão aplicada de 2500 Pa ultrapassa o comprimento utilizado de duto (4 m).

Com a pressão mais baixa (1500 Pa), a partir de 3 m o gradiente de velocidade

não varia mais significativamente ao longo do duto, caracterizando escoamento

desenvolvido. A Figura 48 apresenta os gradientes para os três cenários de pressão

mais baixa e que apresentaram número de Reynolds inferior a 2300 (escoamento

laminar).

Figura 48 Média do gradiente de velocidade axial ao longo do duto

0

5

10

15

20

25

30

35

40

45

50

0

1000

2000

3000

4000

5000

6000

7000

8000

9000

10000

1000 2000 3000 4000 5000 6000 7000

Co

mp

rim

en

to d

e E

ntr

ada

(m)

Nú

me

ro d

e R

eyn

old

s

Pressão de Aplicada (Pa)

Número de Reynolds e Comprimento de Entrada no Regime Permanente

Reynolds

Comprimento

1,E-06

1,E-05

1,E-04

1,E-03

1,E-02

1,E-01

1,E+00

0,00 0,50 1,00 1,50 2,00 2,50 3,00 3,50 4,00

Gra

die

nte

da

velo

cid

ade

(1

/s)

Comprimento (m)

Média do Gradiente de Velocidade no Comprimento

1500Pa

2000Pa

2500Pa

68

5.4 Análises de Sensibilidade

Neste item serão apresentados alguns resultados para diferentes pressões de

entrada e tempos de equilíbrio para ilustrar relações possíveis entre estes

parâmetros e o tempo para o reinício. Apesar das velocidades atingirem valores

elevados a depender da pressão aplicada, o escoamento é laminar pelo menos para

os instantes de tempo iniciais (vide Figura 49) e que são de interesse neste item.

Figura 49 Número de Reynolds medido na saída (x/L=0.975) para as maiores pressões

aplicadas

5.4.1 Sensibilidade à Pressão de Entrada

Diferentes pressões de entrada para uma mesma condição inicial de fluido

em repouso e tempo de equilíbrio (teq igual a 0,6s) foram aplicadas e o tempo

necessário para o reinício foi calculado. As velocidades médias obtidas para os

primeiros instantes são apresentadas na Figura 50.

0

2000

4000

6000

8000

10000

0,00 2,00 4,00 6,00 8,00 10,00 12,00 14,00

Nú

me

ro d

e R

eyn

old

s

Tempo (s)

Número de Reynolds para as Maiores Pressões Medido na Saída (x/L=0.975)

Pressão 6500Pa Pressão 4500Pa Pressão 3500Pa Pressão 2500Pa

69

Figura 50 Velocidade média do fluido em função do tempo e da pressão de entrada para teq

igual a 0,6s

Conforme discutido no item 5.3.1, o comportamento de sólido é

representado por um fluido de viscosidade muito alta. Desta forma, apesar da

Figura indicar que no início o sistema está parado, ao apresentar a mesma

informação em um gráfico log-linear, observa-se que não é este o caso (vide

Figura 51).

Figura 51 Velocidade média do fluido em função do tempo e da pressão de entrada para teq

igual a 0,6s na escala logarítmica

0,00

0,20

0,40

0,60

0,80

1,00

0,00 5,00 10,00 15,00 20,00 25,00

Ve

loci

dad

e M

éd

ia (

m/s

)

Time (s)

Velocidade Média para Diferentes Pressões Aplicadas Medida na Saída (x/L=0.975)

10500Pa

6500Pa

4500Pa

3500Pa

2500Pa

2000Pa

1500Pa

1,0E-06

1,0E-05

1,0E-04

1,0E-03

1,0E-02

1,0E-01

1,0E+00

0,00 5,00 10,00 15,00 20,00 25,00

Ve

loci

dad

e M

éd

ia (

m/s

)

Time (s)

Velocidade Média para Diferentes Pressões Aplicadas Medida na Saída (x/L=0.975) com Escala Log-Linear

10500Pa

6500Pa

4500Pa

3500Pa

2500Pa

2000Pa

1500Pa

70

Como não há um momento no qual a velocidade deixa de ser nula para

caracterizar o início do escoamento, é necessário definir alguma premissa para

calcular o tempo de repartida. O seguinte critério foi utilizado neste e nos demais

itens deste trabalho: calcula-se por interpolação em que tempo a velocidade do

fluido atinge um centésimo de metros por segundo (0,01 m/s). Outra proposta

seria utilizar uma fração da velocidade máxima, entretanto, conforme discutido no

item 5.2.2, pode-se depender de resultados fora do regime laminar.

O resultado obtido para o tempo de repartida em função da pressão aplicada

é apresentado na Figura 52 para o tempo de equilíbrio de 0,6 s e 1,2 s.

Figura 52 Tempo para o reinício em função da pressão de entrada para dois tempos de

equilíbrio

A partir dos resultados da Figura 52, é possível aproximar com boa precisão

o tempo de repartida por uma relação de potência. Observa-se também que quanto

maior a pressão aplicada, menor é o efeito do tempo de equilíbrio e, no limite,

nenhum comportamento tixotrópico seria observado.

5.4.2 Sensibilidade ao Parâmetro b

Os resultados apresentados até o momento consideraram os parâmetros b e

a, da equação 40, iguais a 0,5. A partir da derivada da estrutura do fluido no

tempo para um valor de λeq (estrutura em equilíbrio) qualquer, é possível inferir a

influência de cada um destes parâmetros no desenvolvimento da estrutura.

A Figura 53 ilustra a taxa de quebra do fluido para diferentes valores de b

considerando λeq igual a meio. O parâmetro b afeta tanto a taxa de construção

como a de quebra, mas com maior influência na taxa de quebra.

y = 18951x-1,053

R² = 0,971

y = 47926x-1,077

R² = 0,9802

0,00

5,00

10,00

15,00

20,00

25,00

0 1000 2000 3000 4000 5000 6000 7000

Tem

po

par

a R

ep

arti

da

(s)

Pressão de Entrada (Pa)

Tempo de Repartida por Pressão de Entrada

tEQ = 1.2s

tEQ = 0.6s

71

Figura 53 Derivada do termo tixotrópico (equação 40) para diferentes valores de b

considerando λeq = 0.5

Já o parâmetro a afeta pouco a taxa de quebra ou de construção do fluido. A

Figura 54 ilustra o efeito de a nestas taxas, utilizando a mesma escala que o

gráfico anterior (variando b) para facilitar comparação.

Figura 54 Derivada do termo tixotrópico (equação 40) para diferentes valores de a

considerando λeq = 0.5

Com o objetivo de avaliar qualitativamente a relevância do parâmetro b¸

será apresentada uma breve análise do comportamento de repartida do fluido

mantendo o valor de a igual a 0,5 (caso base). Como as taxas sofrem menor

influência de a do que b, somente a sensibilidade a b será apresentada.

O cenário avaliado considera o sistema partindo do repouso, com fluido

totalmente estruturado, a partir do qual é aplicada uma pressão de 4500 Pa

(cenário de alta pressão, item 5.3.2). Após atingir o regime permanente, a pressão

-6,00

-5,00

-4,00

-3,00

-2,00

-1,00

0,00

1,00

2,00

0,00 0,20 0,40 0,60 0,80 1,00

dλ/

dt

(1/s

)

λ

dλ/dt para diferentes valores de b

b=0.01

b=0.5

b=2

λeq = 0.5

-6,00

-5,00

-4,00

-3,00

-2,00

-1,00

0,00

1,00

2,00

0,00 0,20 0,40 0,60 0,80 1,00

dλ/

dt

(1/s

)

λ

dλ/dt para diferentes valores de a

a=0.01

a=0.5

a=2

λeq = 0.5

72

de entrada é reduzida para 2500 Pa com o objetivo de avaliar a reconstrução do

sistema. Dois valores para b foram avaliados: igual a 0,5 (caso base) e igual a 2.

Quando a pressão de 4500Pa é aplicada, o fluido de maior parâmetro b se

quebra mais rapidamente, e inicia-se o escoamento efetivo antes daquele de menor

valor de b conforme apresentado na Figura 55. No tempo igual a 16s, a pressão de

entrada é reduzida para 2500Pa e a velocidade média é reduzida a uma taxa

aproximadamente igual para ambos os valores de b utilizados. Esse

comportamento é explicado pela proximidade das taxas de construção ilustrada na

Figura 53.

Figura 55 Velocidade média do fluido para a=0,5 e b igual a 0,5 e 2

O comportamento do tempo de reinício em função da pressão para dois

valores de tempo de equilíbrio é apresentada na Figura 56. Observa-se

comportamento semelhante à variação do tempo de equilíbrio (Figura 52) e, por

essa razão, a caracterização do fluido deve ser feito com critério para evitar que

um parâmetro seja ajustado no lugar de outro.

0

1000

2000

3000

4000

5000

0,00

0,50

1,00

1,50

2,00

2,50

0,00 5,00 10,00 15,00 20,00 25,00 30,00

Pre

ssão

de

En

trad

a (P

a)

Ve

loci

dad

e M

éd

ia (

m/s

)

Tempo (s)

Velocidade Média no Tempo Medida na Saída (x/L=0.975) para Alta Pressão de Entrada (4500Pa)

Velocidade b=0.05 Velocidade b=2.00 Pressão Entrada

73

Figura 56 Tempo para o reinício em função da pressão de entrada para dois valores de b

5.4.3 Sensibilidade ao Tempo de Equilíbrio

Além da pressão aplicada, as propriedades tixotrópicas do fluido também

afetam o tempo necessário para o reinício do escoamento. O aumento do tempo de

equilíbrio, por exemplo, apesar de não afetar a velocidade final de escoamento,

reduz a taxa de quebra da estrutura do fluido. Essa redução na velocidade de

quebra implica em maiores tempos de repartida para uma mesma pressão

aplicada. A velocidade para três tempos de equilíbrio é apresentada na Figura 57.

Figura 57 Velocidade média na saída (x/L=0.975) para diferentes tempos de equilíbrio

y = 18951x-1,053

R² = 0,971

y = 16953x-1,124

R² = 0,9673

0,00

2,00

4,00

6,00

8,00

10,00

12,00

0 1000 2000 3000 4000 5000 6000 7000

Tem

po

par

a R

ep

arti

da

(s)

Pressão de Entrada (Pa)

Tempo de Repartida por Pressão para 2 Valores de b

b=2

b=0.5

0,00

0,02

0,04

0,06

0,08

0,10

0,12

0,14

0,16

0,00 20,00 40,00 60,00 80,00

Ve

loci

dad

e M

éd

ia

Tempo (s)

Velocidade Média para Pressão Baixa de Entrada Medida na saída (x/L=0.975)

tEq 1.2s / 1500Pa

tEq 0.6s / 1500Pa

tEq 0.3s / 1500Pa

74

A influência do tempo de equilíbrio (teq) para o tempo de repartida do fluido

é cada vez menor com o aumento da pressão aplicada. Além disso, existe uma

relação linear entre teq e o tempo de repartida, que depende da pressão aplicada,

conforme ilustra a Figura 58. O critério adotado para definir o reinício é o mesmo

do item 5.4.1. Observa-se que a extrapolação das relações obtida para um tempo

de equilíbrio igual a zero é próximo zero, conforme esperado.

Figura 58 Tempo de repartida em função do tempo de equilíbrio para duas pressões de

entrada (baixa e alta)

5.5 Validação Experimental

O objetivo deste item é comparar os resultados apresentados neste trabalho

com os dados experimentais obtidos por Gustavo Moisés (dados não publicados)

para o fluido com os parâmetros reológicos apresentados no item 5.1.

O experimento consistiu em medir a vazão de líquido na saída de um tubo

horizontal de aproximadamente 4m de comprimento e diâmetro interno de

19,05mm a partir de uma pressão de entrada aplicada. A velocidade foi medida

com um velocímetro ultrassônico posicionado próximo a saída do sistema.

Os experimentos utilizados foram aqueles com as pressões na entrada

conforme Figura 59, com as respectivas velocidades obtidas ao longo do tempo

(vide Figura 60).

y = 4,5102x - 0,0949R² = 1

y = 15,636x + 0,485R² = 1

0,00

5,00

10,00

15,00

20,00

25,00

30,00

35,00

0 0,3 0,6 0,9 1,2 1,5 1,8 2,1 2,4

Tem

po

par

a R

ep

arti

da

(s)

Tempo de Equilíbrio (s)

Tempo de Repartida por Tempo de Equilíbrio

P=1500Pa

P=4500Pa

75

Figura 59 Pressão aplicada nos experimentos

Figura 60 Velocidade obtida a partir de dados experimentais

Comportamento semelhante ao experimental foi observado nas simulações

apresentadas, especialmente no item 5.4.1 “Sensibilidade à Pressão de Entrada ”,

onde as velocidades obtidas sem normalização são apresentadas na Figura 61.

0

5000

10000

15000

20000

25000

30000

0,00 50,00 100,00 150,00 200,00

Pre

ssão

(P

a)

Tempo (s)

Resultados Experimentais - Pressão Aplicada

Experimento 1

Experimento 2

Experimento 3

0,00

0,20

0,40

0,60

0,80

1,00

1,20

0,00 50,00 100,00 150,00 200,00

Ve

loci

dad

e (

m/s

)

Tempo (s)

Resultados Experimentais - Velocidade

Experimento 1

Experimento 2

Experimento 3

76

Figura 61 Velocidade média obtida pela simulação numérica

A partir da sensibilidade para o tempo de equilíbrio (item 5.4.3) para baixa

pressão (1500 Pa), com as velocidades normalizadas pela velocidade em regime

permanente, observa-se comportamento semelhante entre a simulação numérica e

o dado experimental, especialmente para tempo de equilíbrio igual a 1,8 segundos

conforme Figura 62.

Figura 62 Comparativo entre dado experimental e simulado

Portanto, as simulações apresentam comportamento coerente com os dados

experimentais utilizados, o que valida qualitativamente o resultado obtido pela

ferramento de CFD e metodologia apresentada.

0,00

0,20

0,40

0,60

0,80

1,00

0,00 5,00 10,00 15,00 20,00 25,00 30,00

Ve

loci

dad

e M

éd

ia (

m/s

)

Time (s)

Velocidade Média para Diferentes Pressões Aplicadas Medida na Saída (x/L=0.975)

2500Pa

2000Pa

1500Pa

-0,20

0,00

0,20

0,40

0,60

0,80

1,00

1,20

0,00 20,00 40,00 60,00 80,00 100,00 120,00

Ve

loci

dad

e /

Ve

loci

dad

e M

áxim

a

Tempo (s)

Comparação entre Dado Experimental e Simulado (medido na saída x/L=0.975)

Experimento 3 Simulação 1500Pa e tEq=0.6s

Simulação 1500Pa e tEq=1.2s Simulação 1500Pa e tEq=1.8s

77

6 Conclusão

O algoritmo e metodologia utilizados para permitir a simulação número de

um fluido tixotrópico na ferramenta de CFD comercial Fluent© foi desenvolvida e

apresentada. Os resultados para o reinício apresentaram comportamento

representativo esperado para fluidos tixotrópicos com semelhança qualitativa com

dados experimentais. Observou-se possível relacionar a pressão aplicada, tempo

de equilíbrio e parâmetro tixotrópico com o tempo necessário para reinício do

escoamento utilizando equações simples, desde que o fluido seja completamente

caracterizado.

Como perspectiva de trabalhos futuros, é possível expandir a simulação e

metodologia apresentadas fluidos compressíveis, adequando o termo da equação

de conservação de quantidade de movimento nos moldes descritos aqui. Além

disso, sugere-se também a avaliação do reinício através do deslocamento por

fluido não tixotrópico o que permitirá melhor compreensão do fenômeno e de

cenários mais alinhados com os esperados na realidade.

Desta forma, com adição desta e outras funcionalidades no algoritmo

apresentado, estudos de engenharia para aplicações reais de repartida de fluido

podem se beneficiar desta ferramenta no futuro.

78

7 Referências bibliográficas

1. CHHABRA, R. P.; RICHARDSON, J. F. Non-Newtonian Flow and Applied

Rheology. Segunda. ed. Oxford, UK: Elsevier, 2008.

2. IRGENS, F. Rheology and Non-Newtonian Fluids. Primeira. ed. [S.l.]:

Springer, 2014.

3. LIMA, D. D. S. Startup Flow of Gelled Crudes after a Shutdown

Comparison between Simulations and Experimental Data. Pontifícia

Universidade Católica do Rio de Janeiro, Departamento de Engenharia

Mecânica. Rio de Janeiro: [s.n.]. 2015. p. 132.

4. BIRD, B. R.; ARMSTRONG, R. C.; HASSAGER, O. Dynamics of Polimeric

Liquids. Segunda. ed. Nova Iorque, NY: John Wiley & Sons, v. I, 1987.

5. SARAMITO, P. A new elastoviscoplastic model based on the Herschel-

Bulkley viscoplastic model. Journal of Non-Newtonian Fluid Mechanics,

Genova, Itália, Maio 2009. 154-161.

6. PAPANASTASIOU, T. C. Flows of Materials with Yield. Journal of

Rheology, v. 31, p. 385-404, 1987. ISSN 10.1122/1.549926.

7. DE SOUZA MENDES, P. R.; DUTRA, E. S. S. Viscosity Function for Yield-

Stress Liquids. Journal of Applied Rheology, v. 14, n. 6, p. 296-302, 2004.

8. ANSYS INC. Ansys Help 15.0. [S.l.]: [s.n.], 2013.

9. DE SOUZA MENDES, P. R. Modeling the thixotropic behavior of structured

fluids. Journal of Non-Newtonian Fluid Mechanics, v. 164, p. 66-75, 2009.

10. FONSECA, C. et al. Flow of Elasto-Viscoplastic Thixotropic Liquids Past a

Confined Cylinder. Journal of Non-Newtonian Fluid Mechanics, v. 193, p.

80-88, 2013.

11. GOLDBERG, D. What Every Computer Scientist Should Know About

Floating-Point Arithmetic. Computing Surveys, 1991.

12. BARNES, H. A. Thixotropy - a review. Journal of Non-Newtonian Fluid

Mechanics, v. 70, p. 1-33, Janeiro 1997.

13. AHMADPOUR, A. Restart Problem of Waxy Crude Oil: a Numerical

79

Study. [S.l.]: University of Tehran, 2014.

14. BAUER, W. H.; COLLINS, E. A. Thixotropy and Dilatancy. In: EIRICH, F.

R. Rheology Theory and Applications. Nova york, NY: Academic Press, v.

4, 1967. Cap. 8.

15. INTERNATIONAL UNION OF PURE AND APPLIED CHEMISTRY. Work

Softening. IUPAC Compendium of Chemical Terminology - the Gold

Book, p. 1217, 2014. Disponivel em:

<http://goldbook.iupac.org/W06691.html>. Acesso em: 25 Julho 2015.

16. GREEN, H.; WELTMANN, R. N. Analysis of the thixotropy of pigment-

vehicle suspesions. Industrial and Engineering Chemistry, Analytical

Edition, Easton, Pa., v. 15, p. 201-206, 1943. ISSN 0096-4484.

17. MEWIS, J.; WAGNER, N. J. Thixotropy. Advances in Colloid and Iterface

Science, v. 147-148, p. 214-227, March-June 2009.

18. MUJUMDAR, A.; BERIS, A. N.; METZNER, A. B. Transient Phenomena in

Thixotropic Systes. Journal of Non-Newtonian Fluid Mechanics, v. 102, p.

157-178, 2002. ISSN S0377-0257(01)00176-8.

19. POTANIN, A. A. et al. Microrheological Modeling of Weakly Aggregated

Dispersions. The Journal of Chemicall Physics, v. 4, n. 102, p. 5845-5853,

1995. ISSN DOI: 10.1063/1.469317.

20. MOORE, F. The Rheology of Ceramic Slips and Bodies. Transactions of the

British Ceramic Society, v. 58, p. 470-494, 1959.

21. LIVESCU, S. Mathematical modeling of thixotropic drilling mud and crude

oil flow in wells and pipelines - A Review. Journal of Petroleum Science

and Engineering, v. 98-99, p. 174-184, 2012.

22. HOUSKA, M. Engineering aspects of the rheology of thixotropic liquids.

PhD Thesis. ed. [S.l.]: Czech Technical University ofPrague, 1981.

23. WACHS, A.; VINAY, G.; FRIGAARD, I. A 1.5D numerical model for the

start up of weakly compressible flow of a viscoplastic and thixotropic fluid in

pipelines. Journal of Non-Newtonian Fluid Mechanics, v. 159, p. 81-94,

2009.

24. AHMADPOUR, A.; SADEGHY, K. Start-up flows of Dullaert-Mewis

viscoplastic-thixoelastic fluids: a two dimensional analysis. Journal of Non-

80

Newtonian Fluid Mechanics, v. 214, p. 1-17, 2014.

25. DE SOUZA MENDES, P. R. Thixotropic elasto-viscoplastic model for

structured fluids. Soft Matter, n. 7, p. 2471–2483, 2011. ISSN DOI:

10.1039/c0sm01021a.

26. DE SOUZA MENDES, P. R. et al. Startup flow of gelled crudes in pipelines.

Journal of Non-Newtonian Fluid Mechanics, n. 179-180, p. 23-31, 2012.

ISSN DOI: 10.1016/j.jnnfm.2012.05.003.

27. DULLAERT, K.; MEWIS, J. A scrutctural kinetics model for thixotropy.

Journal of Non-Newtonian Fluid Mechanics, n. 139, p. 21-30, 2006.

28. POTANIN, A. 3D simulations of the flow of thixotropic fluids, in large-gap

Couette and vane-cup geometries. Journal of Non-Newtonian Fluid

Mechanics, v. 165, p. 299-312, 2010.

29. AZEVEDO, L. F. Mecânica dos fluidos I - Notas de Aula. 02 de Agosto de

2010. ed. Rio de Janeiro: Pontifícia Universidade Católica do Rio de Janeiro,

Departamento de Engenharia Mecânica, 2010.

30. PATANKAR, S. V. Numerical Heat Transfer and Fluid Flow. Nova Iroque,

NY: Hemisphere Publishing, 1980.

31. BAKKER, A. Lecture 5 - Solution Methods. Computational Fluid

Dynamics, 2006. Disponivel em:

<http://www.bakker.org/dartmouth06/engs150/05-solv.pdf>. Acesso em: 01

jun. 2015.

32. DE SOUZA MENDES, P. R.; THOMPSON, R. L. A Unified Approach to

Model Elasto-viscoplastic Thixotropic Yield-stress Materials and Apparent

Yield-stress Fluids. Rheol Acta, Berlin, 23 Abril 2013. 22.

33. LE MEINS, J. F.; MOLDENAERS, P.; J., M. Suspensions in polymer melts.

Industrial & Engineering Chemistry Research, v. 41, p. 6297-6304, 2002.

34. GOODEVE, C. F. A general theory of thixotropy and viscosity. Transactions

of the Faraday Society, v. 35, p. 342-358, 1939.