2.1Sobre as modelagens deste relatório e outras

Neste relatório, descreve-se o desenvolvimento de modelagens computacionais através do SisBaHiA^® - Sistema Base de Hidrodinâmica Ambiental, em um estudo da circulação hidrodinâmica no Estuário do Rio Paraíba do Norte, PB, incluindo análise de tempos hidráulicos característicos, como tempos de taxas de renovação e idade da água. Embora apenas casos simples de Modelos Hidrodinâmico e de Transporte Euleriano do SisBaHiA^® estejam exemplificados, o interessado pode desenvolver outros modelos a partir da base de dados disponibilizada, como descrito na Ref. Técnica do SisBaHiA^® , v. ROSMAN (2018):

• Modelos de Transporte de Sedimentos Coesivos e Não Coesivos e Evolução Morfológica do Fundo;
• Modelos de Qualidade de Água e Eutrofização;
• Modelos de Transporte Lagrangeanos Determinístico e Probabilístico;
• Modelo de Campo Próximo e Decaimento Bacteriano para Emissários e Fontes Pontuais
• Modelo de Geração de Ondas;
• Modelo de Propagação de Ondas.

Dependendo dos efeitos incluídos na modelagem hidrodinâmica, outros modelos podem ser acoplados. Seguimos a reproduzir a introdução do Capítulo 3 da Ref. Técnica do SisBaHiA^® , v. ROSMAN (2018), sobre a aplicação e acoplamento dos modelos:

• MH = Modelos Hidrodinâmicos com ou sem efeitos de ondas de curto período:
• MGO = modelo de geração de ondas por vento local no domínio de modelagem;
• MPO = modelo de propagação de ondas remotas geradas fora do domínio de modelagem do MH.
• MQA = Modelos de Qualidade de Água, v. seção 6 da Ref. Técnica do SisBaHiA^® (ROSMAN, 2018), para transporte de constituintes que, dependendo do caso, têm comportamento ativo: salinidade e calor (temperatura). Os demais constituintes, e.g. OD, DBO e nutrientes, sempre têm comportamento passivo e por isso são usualmente computados de modo desacoplado.
• MS = Modelos Sedimentológicos, veja capítulo 4 da Ref. Técnica do SisBaHiA^® (ROSMAN, 2018) para sedimentos não-coesivos e capítulo 5 da Ref. Técnica do SisBaHiA^® (ROSMAN, 2018) para sedimentos coesivos, em especial seção 5.2.

As seguintes opções são possíveis, envolvendo acoplamento ou não entre os modelos. Quando modelos rodam acoplados são interdependentes. Note que em todos os casos abaixo MH significa modelo hidrodinâmico com ou sem efeitos de ondas de curto período (MGO e/ou MPO):

1. MH: Análises de circulação hidrodinâmica sem efeitos de gradientes de densidade e sem evolução morfodinâmica de batimetria. Neste caso o MH definido roda apenas para gerar padrões de circulação hidrodinâmica com ou sem efeitos de ondas de curto período. A inclusão de efeitos de ondas de curto período pode ser feita por acoplamento de um MGO e/ou de um MPO. A inclusão de efeitos de onda no MH altera a circulação hidrodinâmica, principalmente em áreas mais rasas e perto da linha de costa, pois afeta as tensões de atrito no fundo e inclui efeitos de tensões de radiação geradas por ondas. Tais efeitos são especialmente relevantes para processos sedimentológicos em zonas costeiras com praias, em lagos, em baías, e similares, sendo irrelevantes em rios. Cenários de MQA e MS podem ser vinculados ao MH desenvolvido, mas nesta opção rodariam desacoplados.
2. MH+MQA: Análises de circulação hidrodinâmica com efeitos de gradientes de densidade e sem evolução morfodinâmica de batimetria. Neste teremos MH+MQA acoplados, o MH definido, com ou sem efeitos de ondas, roda junto com um MQA acoplado para transporte de sal e/ou temperatura. Outros cenários de MQA e MS podem ser vinculados ao MH desenvolvido, mas nesta opção rodariam desacoplados.
3. MH+MS: Análises de circulação hidrodinâmica sem efeitos de gradientes de densidade e com evolução morfodinâmica de batimetria. Neste caso teremos MH+MS acoplados, o MH desenvolvido, com ou sem efeitos de ondas, roda junto com MS de sedimentos não-coesivos e/ou coesivos acoplados para transporte de sedimentos e processos de erosão e assoreamento. Outros cenários de MQA e MS podem ser vinculados ao MH desenvolvido, mas nesta opção rodariam desacoplados.
4. MH+MQA+MS: Análises de circulação hidrodinâmica com efeitos de gradientes de densidade e com efeitos de evolução morfodinâmica. Neste caso teremos MH + MQA + MS acoplados, o MH definido, com ou sem efeitos de ondas, roda junto com um MQA acoplado para transporte de sal e/ou temperatura e junto com MS de sedimentos não-coesivos e/ou coesivos acoplados para transporte de sedimentos e processos de erosão e assoreamento. Outros cenários de MQA e MS podem ser vinculados ao MH desenvolvido neste caso, mas rodariam desacoplados.

Vale destacar que no SisBaHiA^® pode-se rodar de modo acoplado MH + MGO + MPO + MQA(Salinidade) + MQA(Temperatura) + MS(não-coesivos) + MS(coesivos). Neste caso, os modelos sentem efeitos uns dos outros, pois são interdependentes.

Em qualquer dos casos, outras análises envolvendo aspectos de qualidade de água ou transporte de sedimentos a serem modeladas com um MQA ou um MS podem ser feitas de modo desacoplado. Nestes casos, os MQA e MS serão modelos clientes vinculados a um MH e só podem rodar depois de obtidos os resultados do MH vinculante.

O modelo hidrodinâmico do SisBaHiA^® é chamado de FIST3D (filtered in space and time 3D)Rosman, P.C.C. – Modeling Shallow Water Bodies via Filtering Techniques. Ph.D. thesis Civil Engineering Department, Massachusetts Institute of Technology, 1987 . Trata-se de um eficiente modelo numérico hidrodinâmico, tridimensional, para escoamentos homogêneos e de grande escala. O modelo é composto por dois módulos:

1. Um módulo 2DH, que calcula correntes promediadas na vertical e a elevação da superfície livre. Por razões de eficiência numérica, o modelo FIST3D sempre inclui esse módulo.
2. Um módulo 3D que calcula o campo de velocidades tridimensional, com duas opções de metodologia definidas pelo usuário.

Em modelagens hidrodinâmicas 3D, os dois módulos funcionam acoplados e são interdependentes. O módulo 2DH pode ser executado como um módulo independente, se o modelador desejar apenas valores promediados na vertical. Entretanto, como já mencionado, se um campo de escoamento 3D for desejado, existem duas opções disponíveis, e ambas são acopladas ao módulo 2DH:

• a) Modelo 3D completo, totalmente numérico
• b) Modelo 3D analítico-numérico para perfis de velocidade no campo de escoamento horizontal. Esta opção é mais eficiente em termos computacionais, mas apenas inclui a aceleração advectiva do módulo 2DH. Portanto, os resultados obtidos são menos precisos em regiões nas quais as acelerações advectivas variem significativamente ao longo da profundidade. Nessa opção, os perfis de velocidade são computados através de uma solução que é função das velocidades 2DH promediadas na vertical, elevação da superfície livre, rugosidade equivalente de fundo do módulo 2DH, e da velocidade do vento atuando na superfície livre da água.

Em ambos os casos os gradientes de elevação da superfície e as tensões de atrito no fundo acoplam os dois módulos, sendo que um acoplamento adicional é garantido ao se forçar que as velocidades médias na vertical sejam iguais em ambos os módulos. O módulo 3D utiliza as elevações da superfície livre obtidas no módulo 2DH. Em sequência, o módulo 2DH utiliza as tensões de atrito no fundo obtidas dos perfis verticais de velocidade calculadas no módulo 3D. No caso da opção estritamente numérica, (a), os perfis de velocidades calculados no módulo 3D são integrados na vertical e comparados com as velocidades 2DH obtidas no módulo 2DH. Caso necessário, os perfis de velocidade 3D são ajustados de modo a fazer coincidir as velocidades médias na vertical nos dois módulos.

4.1Definição do domínio de interesse

O domínio de modelagem é limitado por contornos terrestres e abertos. Na definição do domínio, em geral, os contornos abertos são limites de conveniência entre massas de água. No caso do ERPN exemplificado, ambos os contornos foram definidos a partir de imagens do Google Earth e foram unidos e redesenhados usando o programa Surfer, constituindo o mapa base para este estudo, cf. Figura 2. Na figura, o mapa ficou propositalmente maior que o domínio, de modo a facilitar a localização e ter espaço para incluir escalas e outras informações usuais em mapas para relatórios técnicos.

O contorno de terra representa as margens de baías ou lagoas, seções de rios ou canais e margens de ilhas que estão dentro do domínio de modelagem. Já o contorno de água delimita o domínio de modelagem não prescrito pelo contorno físico, como é o caso de entradas de baías e estuários. O alcance de ambos contornos devem ser pensados de acordo com o interesse específico de cada modelagem.

4.2Malha de discretização

Após exportar o mapa base e seus respectivos contornos para o SisBaHiA^®, foi confeccionada a malha de elementos finitos utilizada na discretização do domínio espacial do modelo, Figura 3. A malha do domínio resultou em 2302 elementos totais, sendo 1888 elementos quadrangulares e 127 triangulares, e contendo um total de 10504 nós de cálculo, como mostra o quadro de informações na tela inicial do item “Malhas & Domínios de Modelagem”, Figura 4. Esta tela reúne estas informações assim que a malha é criada no item “Modelagem” encontrado no canto superior esquerdo da tela inicial do SisBaHiA^®.

4.3Batimetria

Os dados de batimetria do ERPN adotados foram provenientes de cartas náuticas da Diretoria de Hidrografia e Navegação da Marinha do Brasil (DHN). Especificamente, das cartas náuticas Porto de Cabedelo no 830 (escala 1:15.000) e Proximidades do Porto de Cabedelo nº 806 (escala 1:50.000).

Os dados de batimetria foram reunidos e organizados para serem posteriormente interpolados através do SisBaHiA^®, que utiliza o programa Surfer. Neste caso, adotou-se o método de interpolação Kriging. A Figura 5 apresenta o mapa com as batimetrias resultantes para o domínio de modelagem.

Os dados interpolados de batimetria utilizados na modelagem podem ser extraídos acessando a base de dados no SisBaHiA^® em Malha & Domínios de Modelagem na aba Nós, e clicando posteriormente com o botão direito na tabela mostrada na Figura 6, escolhendo copiar tabela.

4.4Rugosidade equivalente de fundo

A tensão de atrito no fundo depende da amplitude da rugosidade equivalente de fundo, que é determinada em função da natureza do sedimento da área a ser modelada. Os valores desta rugosidade devem ser dados ao modelo hidrodinâmico em todo seu domínio. Estes são adotados através de uma aproximação seguindo os valores sugeridos por ABBOTT e BASCO (1989) adaptado por ROSMAN (2018), como mostrado na Tabela 1.

Os valores adotados para a amplitude da rugosidade equivalente de fundo, sem os efeitos de ondas, foram inferidos a partir da distribuição granulométrica do fundo. A Figura 7 mostra a configuração da rugosidade do fundo adotada para o modelo hidrodinâmico.

Resultados obtidos em modelagem hidrodinâmica preliminar, bem como informações disponíveis na carta náutica 830, sugerem que no sul do estuário predomina a presença de siltes, enquanto na região central, até sua embocadura, há predominância de areias médias. Nas praias adjacentes à saída do estuário encontram-se fundos compostos por areias médias e grossas, além da presença de recifes de corais na região.

Os dados interpolados de rugosidade equivalente de fundo utilizados na modelagem podem ser extraídos acessando a base de dados no SisBaHiA^® em Malha & Domínios de Modelagem na aba Nós, e clicando posteriormente com o botão direito na tabela mostrada na Figura 6

5.1Cenários de Modelagem

Para uma melhor caracterização da hidrodinâmica e da renovação das águas do estuário, é interessante apresentar análises com diferenças sazonais que consideram diferentes padrões de ventos, marés e descargas fluviais. Sendo assim, foram propostos dois cenários de simulação. O primeiro cenário proposto foi o representativo de cheia, incluindo os meses de março, abril e maio de 2017, o segundo cenário adotado foi o representativo de seca, incluindo neste, os meses de outubro, novembro e dezembro de 2017. A escolha de tais períodos de simulação ocorreu pelo fato da região não ter uma variabilidade climática distinta entre as estações do ano, assim, a maior representatividade dos períodos definidos se deu de acordo com os meses de maiores e menores vazões dos afluentes. O ano escolhido está de acordo com a disponibilidade de dados mais recentes.

A Figura 8 mostra a tela inicial do modelo hidrodinâmico do SisBaHiA^® com as configurações iniciais da simulação de cheia, sendo as mesmas aplicadas para o modelo de seca. O passo a passo para a montagem de um modelo está presente nos próximos itens de forma simplificada. Para maiores detalhes desta montagem e também especificações de ferramentas presentes neste modelo acesse o link: http://www.sisbahia.coppe.ufrj.br/ManualSisbahia.htm, do “Manual do Usuário do SisBaHiA^®.

5.2Dados oceanográficos

Os itens a seguir descrevem os dados oceanográficos inseridos no modelo bem como explica a metodologia utilizada para a melhor representação da hidrodinâmica no domínio de modelagem.

5.2.1Marés

Os níveis de maré inseridos como condicionantes na fronteira aberta com o mar foram obtidos a partir do conjunto de constantes harmônicas do Porto de Cabedelo, disponibilizadas pela Fundação FEMAR e expostas em ordem de magnitude na Tabela 2. Tais constantes são referidas a medições de maré durante um período de 365 dias entre 31 de janeiro de 1981 e 31 de janeiro de 1982. Ressalta-se que a previsão de maré obtida a partir das constantes harmônicas da FEMAR foi comparada com uma série gerada a partir de constantes harmônicas do modelo global FES2014 (Finite Element Solution). Os resultados foram praticamente idênticos. Dado o resultado obtido por meio desta comparação, bem como o período de 365 dias de medição, assumiu-se que o dado de maré considerando as constantes harmônicas da FEMAR é confiável.

A maré na região modelada é semidiurna com desigualdade diurna, isto é, ocorrem diferenças de altura entre duas baixa-mares ou entre duas preamares no mesmo dia. Durante a sizígia, as alturas de maré são frequentemente maiores que 2,5 m, conforme exposto na Figura 9.

Para extrair os dados de elevação fornecidos para a fronteira aberta nesta modelagem, basta abrir a base de dados no SisBaHiA, acessar o Modelo Hidrodinâmico do período de Cheia ou de Seca, na aba Parâmetros, e clicar com o botão esquerdo no botão com a lupa mostrado na Figura 10 . A tabela de constantes harmônicas está disponível na aba “Constantes Harmônicas” do menu Malhas e Domínios de Modelagem.

5.3Dados hidro-meteorológicos

Este item reúne dados sobre vazões e ventos inseridos para a modelagem hidrodinâmica e suas respectivas fontes e metodologia aplicada.

5.3.1Vazões

A vazão fluvial mais significativa que aflui ao ERPN é aquela proveniente desse próprio rio. O modelo considerou dados de vazão fluvial média mensal de todos os afluentes. As vazões médias de cada mês do Rio Paraíba foram adaptadas de XAVIER et al. (2012), que compilou dados da série histórica do HIDROWEB, sistema de informações hidrológicas da Agência Nacional de Águas (ANA, 2016), para a bacia hidrográfica do Rio Paraíba, considerando o período entre 1970 e 2006.

A estação fluviométrica instalada mais próxima da região de interesse é a estação Ponte da Batalha, responsável por drenar uma bacia hidrográfica com área de 19.000 km², segundo informação da ANA. No entanto, os dados disponíveis para tal estação apresentavam valores inconsistentes. À montante desta estação, têm-se a estação Guarita, que drena uma bacia hidrográfica com área de 17.400 km², valor que é menos de 10% menor que a área drenada pela estação Ponte da Batalha. Assim sendo, a vazão média mensal dada na fronteira fluvial do Rio Paraíba foi considerada como 10% maior que os dados calculados para a estação Guarita.

Devido à escassez de dados, por não haver estações fluviométricas suficientes e medições na região, as vazões fluviais dos demais rios contribuintes ao ERPN foram inferidas e inseridas na modelagem utilizando valores de fluxo médio mensal baseados e comparados com as vazões do Rio Paraíba do Norte, conforme Tabela 3. Os rios considerados na modelagem estão desenhados no mapa base da Figura 2.

Para introduzir os rios na modelagem utilizou-se o tipo de nó “V9” que permite modelar, através de um talude prescrito, o efeito da maré sobre o fluxo do rio, invertendo o sentido do fluxo ou diminuindo a sua velocidade. A aba no qual está informação é definida está disposta na imagem da tela do SisBaHiA^® da Figura 12. A intrusão deste efeito da maré em rios é importante quando o efeito da maré for além do ponto final do rio a montante considerado na modelagem, como ocorreu no caso que está sendo apresentado neste relatório. Este tipo de nó fornece dados variáveis ao longo do tempo, garantindo o efeito de sazonalidade do estudo.

Neste trabalho foram definidos taludes de forma preliminar. Na Ref. Técnica do SisBaHiA^® essa condição é explicada em maior detalhe, sendo apresentada metodologia de estimativa de taludes.

No SisBaHiA^®, dados de vazão devem ser prescritos em m³/s/m, distribuídos em nós de um segmento de fronteira. A vazão é no modelo como uma distribuição parabólica quadrática pelo lado do respectivo elemento, centrada no meio do lado do elemento. Os rios afluentes ao ERPN foram representados com apenas um elemento na seção de montante. A vazão nodal foi distribuída de forma parabólica, com valores zero nas margens e q = 1.5Q/L no ponto central, onde q é a vazão distribuída no lado do elemento (m³/s/m), Q é a vazão fluvial (m³/s) e L é a largura da seção transversal do elemento (m). Essas informações são definidas na aba Fronteiras/Contornos/Terrestres, como mostra a imagem da tela do SisBaHiA^® na Figura 12.

Como percebe-se na Figura 12, que mostra os valores de vazão por metro linear prescritos para o nó central da seção de montante do Rio Paraíba do Norte, os valores de vazão são acompanhados por um sinal negativo, que indica que o fluxo está entrando no domínio, ou seja, é um afluxo. Se fosse positivo, indicaria que o fluxo estaria saindo do sistema, caracterizando um efluxo. Esse sinal é sempre referente ao sentido do vetor normal à fronteira, que aponta para o exterior da malha, no respectivo lado. Maiores explicações podem ser encontradas na Ref. Técnica do SisBaHiA^®, v. Rosman (2018).

5.3.2Ventos

A tensão de atrito do vento na superfície livre é calculada pelo modelo quando se prescreve valores de velocidades em estações localizadas sobre o domínio de modelagem, e o modelo interpola valores de direção e intensidade dos ventos para cada nó da malha. Com o objetivo de introduzir na modelagem um vento variável no tempo e variado no espaço, adotou-se séries temporais de seis estações, como mostra a tela do SisBaHiA^® na Figura 13. As estações, ilustradas na Figura 14, possuem séries temporais de intensidade de velocidade e direção para o tempo total de simulação.

Os dados de ventos utilizados neste trabalho são dados de reanálise disponibilizados pelo site http://apps.ecmwf.int/datasets/data/interim-full-daily/levtype=sfc/, do ECMWF – European Comunity Medium scale Weather Forecast. Os dados utilizados, obtidos em grade regular de 0.125° x 0,125° com resolução temporal de 3h, representam o vento a 10 m de altura. As séries temporais obtidas contem valores a cada seis horas de velocidades das componentes Leste (U) e Norte(V) dos ventos. As séries temporais de vento obtidas em cada estação foram interpoladas para os nós da malha do modelo hidrodinâmico, a fim gerar o campo de vento para todo o domínio.

Repare na Figura 13, que os dados inseridos no modelo estão na forma de magnitude de velocidade e direção do vento em ângulo azimute.

6Modelo de transporte Euleriano

O cálculo da taxa de renovação e idade da água foi realizado através do Modelo de Transporte Euleriano do SisBaHiA^®. Este pode ser utilizado para análises gerais de qualidade de água alterada por algum escalar cuja concentração não interfere na hidrodinâmica. O modelo admite ampla variedade de tipos de reações cinéticas entre a substância e água do meio receptor.

Foram modelados dois cenários para cada uma das análises. Um cenário de cheia, calculado a partir da modelagem hidrodinâmica do meses de março, abril e maio de 2017, e outro cenário de seca, calculado a partir de outubro, novembro e dezembro de 2017. A imagem da Figura 20, mostra a tela inicial deste modelo com as especificações da modelagem representativa de cheia da taxa de renovação da água.

6.1Taxa de renovação

Para preparar o modelo da taxa de renovação marca-se no instante inicial um valor de referência igual a 0% para as águas que estão dentro do sistema. A renovação devida ao efeito da maré considera que a água que entra no domínio pela fronteira aberta possui valor de referência igual a 100%. Para o aporte fluvial, considera-se que as vazões fluviais afluentes ao domínio, na cabeceira dos rios, possuem o valor 100% também. A renovação total das águas é calculada a partir da combinação destas duas condições iniciais.

Desta forma, a porcentagem de água renovada no domínio será diretamente proporcional a esta concentração, variando de 0% a 100%, sofrendo apenas advecção e dispersão em todo domínio. Ao longo do tempo de simulação, trocas de águas e processos de mistura ocorrem desde o instante inicial, originadas das bacias hidrográficas e da região costeira. Os valores resultantes destas simulações são interpretados como porcentagem de água nova que entrou na região de interesse.

É válido ressaltar que as águas consideradas “novas” não são necessariamente águas limpas. O modelo em questão não representa a qualidade das águas, pois não foram prescritas quantidades de matéria orgânica ou nutrientes afluentes no sistema, por exemplo. Mais explicações sobre a concepção do modelo de Taxa de Renovação podem ser encontradas no capítulo 5 da Ref. Técnica do SisBaHiA, v. Rosman (2018).

A Figura 21 mostra a condição inicial dada ao Modelo de Transporte Euleriano que representa a taxa de renovação espacial do domínio de modelagem. Foram rodados cenários de cheia e seca durante noventa dias, que correspondem aos dias de simulação do modelo hidrodinâmico.

6.2Idade da água

A análise da idade da água permite analisar quanto tempo a água permanece em diferentes setores ao longo do tempo. Este tempo é estimado a partir do decaimento de uma substância passiva marcadora de idade presente na água. Para poder determinar o tempo de decaimento, é obrigatório que a substância marcadora de idade tenha reação cinética de decaimento de primeira ordem, com taxa constante de k > 0, sem outros efeitos de perdas e ganhos de massa.

Para conceituar o cálculo, é considerado um volume de água bem misturado com concentração inicial C₀ de tal substância. Sendo a cinética de decaimento de primeira ordem, a variação no tempo da concentração da substância marcadora de idade, C_(t)0, é dada por: dC⁄dt = -kC, cuja solução analítica leva a: C(t) = C₀ exp⁡(-kt) e, portanto, t = - ln⁡ (C⁄C₀).

Conhecida a concentração inicial C₀ e uma concentração C registrada posteriormente, determina-se diretamente o tempo de decaimento decorrido entre o instante inicial e o instante do registro de C. A diferença entre estes dois instantes define a “Idade da Água” no instante do registro.

Como se trata de um corpo de água natural, com entradas de água nova com C = C₀ em diferentes pontos e tempos, bem como saídas de água que já sofreram decaimento e que apresentam concentrações diferentes em diversos pontos e tempos, a concentração é função do tempo e do espaço, C (x,y,t). Consequentemente, a Idade da Água será uma função variando no espaço e no tempo: I_A (x,y,t) = ( -ln ⁡(C (x,y,t)) / C₀) ⁄ k.

Pelo exposto, para preparar o modelo de idade da água marca-se no instante inicial C(x,y,t₀) = C₀ = 1 em todo o domínio de modelagem, Figura 22. Consequentemente, no início I_A = 0, pois o ln (1) = 0 em todos os locais. As novas águas que entram no domínio pelos principais afluentes tem idade zero e, por isso, devem ter concentração de substância marcadora de idade igual a um.

A medida em que as águas iniciais e as águas novas com I_A = 0 vão se misturando no domínio de interesse e sendo transportadas, o valor de C vai diminuindo em cada local em função do processo de decaimento. Com isso, o valor de I_A torna-se diferente em cada ponto, pois depende da magnitude das correntes e da turbulência no local.