Joseph Harari
fundamentos de
MODELAGEM NUMÉRICA
em OCEANOGRAFIA
2015
São Paulo
fundamentos de
MODELAGEM NUMÉRICA
em OCEANOGRAFIA
Joseph Harari
2015
Editoração
SALT | Sea & Limno Technology
Diagramação
Danilo Rodrigues Vieira
Thiago Marques Coelho
Capa:
Leandro Inoe Coelho
Catalogação na Publicação (CIP)
Ficha Catalográfica feita pelo autor
H254f Harari, Joseph
Fundamentos de modelagem numérica em
Oceanografia / Joseph Harari. – São Paulo, 2015
246 p.; il., color; 21, 59 × 27, 94 cm
ISBN 978-85-918934-0-9
1. Modelagem numérica. 2. Oceanografia.
I. Tı́tulo.
CDD: 551.46
CDU: 556.5
Sumário
Sumário
i
Lista de Figuras
vi
Lista de Tabelas
xii
Prefácio
xiii
Introdução
xiv
1 Conceitos básicos em Modelagem Numérica
1.1 Introdução . . . . . . . . . . . . . . . . . . .
1.2 Método de diferenças finitas . . . . . . . . .
1.3 Diferenças finitas de alta ordem . . . . . . .
1.4 Equação da advecção . . . . . . . . . . . . .
1.5 Estrutura de um modelo . . . . . . . . . . .
1.6 Método de elementos finitos . . . . . . . . .
1.7 Validação dos resultados de um modelo . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
1
1
4
7
8
9
10
13
2 Esquemas de diferenças finitas explı́citos, implı́citos e iterativos
2.1 Equação da difusão . . . . . . . . . . . . . . . . . . . . . . .
2.2 Esquemas de diferenças finitas explı́citos, implı́citos e iterativos
2.3 Solução de esquemas implı́citos . . . . . . . . . . . . . . . .
2.4 Complementação . . . . . . . . . . . . . . . . . . . . . . . .
17
17
18
20
22
3 Análise de erros em modelos numéricos
25
i
Sumário
3.1
3.2
3.3
3.4
3.5
3.6
3.7
3.8
Conceitos básicos . . . . . . . . . . . . . . . . . . . . . . . .
Método de von Neumann para análise de estabilidade . . . .
Análise de estabilidade para um esquema com 3 nı́veis de tempo
Análise de estabilidade para um esquema implı́cito . . . . .
Erros de dispersão computacional . . . . . . . . . . . . . . .
Erros associados ao modo computacional (ruı́do) . . . . . . .
Erros associados à difusão numérica . . . . . . . . . . . . . .
Conclusão . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4 Equações e suas condições de estabilidade. Equações combinadas e formulações 2D e 3D
4.1 Resumo das condições de estabilidade para soluções da equação
da advecção . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2 Resumo das condições de estabilidade para soluções da equação
da difusão . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.3 Equação do decaimento—o método direto de análise de estabilidade . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.4 Equação do oscilador . . . . . . . . . . . . . . . . . . . . . .
4.5 Difusão e decaimento . . . . . . . . . . . . . . . . . . . . . .
4.6 Soluções para a equação da advecção–difusão–decaimento 1D
4.7 Exemplos de discretização de equações 2D e 3D . . . . . . .
25
26
30
31
33
35
35
36
44
44
44
45
47
48
49
50
5 Condições de contorno computacionais
5.1 Motivação básica . . . . . . . . . . . . . . . . . . . . . . . .
5.2 Possı́veis soluções computacionais para os contornos . . . . .
5.3 Consequências das soluções adotadas . . . . . . . . . . . . .
5.4 Modelagem considerando contornos continentais e ilhas . . .
5.5 Modelagem da dispersão para uma região geográfica especı́fica
58
58
59
62
66
68
6 Filtragens e instabilidade não linear
6.1 Filtragem de resultados de modelo (no espaço) .
6.2 Filtragem no tempo . . . . . . . . . . . . . . . .
6.3 Instabilidade não-linear . . . . . . . . . . . . . .
6.4 Possı́veis soluções para a instabilidade não-linear
73
73
75
77
80
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
7 Esquemas numéricos de alta precisão
87
7.1 Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87
7.2 Equação da advecção . . . . . . . . . . . . . . . . . . . . . . 87
ii
Sumário
7.3
7.4
Equação da difusão . . . . . . . . . . . . . . . . . . . . . . . 98
“Splitting” (divisão ou parcelamento) . . . . . . . . . . . . . 105
8 Equações Elı́pticas
8.1 Introdução . . . . . . . . . . . .
8.2 Método da relaxação . . . . . .
8.3 Relaxação sequencial . . . . . .
8.4 Super-relaxação . . . . . . . . .
8.5 Métodos diretos . . . . . . . . .
8.6 Método da Eliminação de Gauss
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
107
107
108
109
110
110
111
9 Modelos numéricos hidrodinâmicos 1D
113
9.1 Solução do sistema de equações de águas rasas 1D simplificado
(com condições de contorno) . . . . . . . . . . . . . . . . . . 113
9.2 Solução explı́cita de primeira ordem para o sistema hidrodinâmico 1D . . . . . . . . . . . . . . . . . . . . . . . . . . . 117
9.3 Solução explı́cita de segunda ordem para o sistema hidrodinâmico 1D . . . . . . . . . . . . . . . . . . . . . . . . . . . 120
9.4 Soluções implı́cita e semi-implı́cita para o sistema hidrodinâmico 1D . . . . . . . . . . . . . . . . . . . . . . . . . . . 122
9.5 Solução do sistema de equações de águas rasas 1D com dissipação125
10 Modelos Numéricos Hidrodinâmicos 2D
128
10.1 Equações hidrodinâmicas básicas 2D . . . . . . . . . . . . . 128
10.2 Modelo hidrodinâmico 2D com solução de primeira ordem no
tempo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 129
10.3 Modelo hidrodinâmico 2D com solução de segunda ordem no
tempo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132
10.4 Soluções implı́citas e semi-implı́citas para os termos de decaimento e difusão (com opção não linear para o decaimento) . 133
10.5 Imposição de condições iniciais e de contorno e a energia do
sistema . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 135
10.6 Condições de estabilidade para a modelagem hidrodinâmica 2D137
10.7 Análise da dispersão para a modelagem 2D . . . . . . . . . . 140
10.8 Modelagem hidrodinâmica 2D considerando contornos continentais e ilhas . . . . . . . . . . . . . . . . . . . . . . . . . . 142
10.9 Equações básicas de modelos 2D simplificados . . . . . . . . 142
10.10Modelagem 2D com a equação da vorticidade . . . . . . . . 144
iii
Sumário
10.11Discretização da modelagem 2D com a equação da vorticidade148
11 Gradeamento em modelos 2D
11.1 Grades alternadas (de Arakawa) . . . . . .
11.2 Refinamento e aninhamento de grades . . .
11.3 Tipos de grades: inclinadas, curvilı́neas,
contorno terrestre variável no tempo . . .
150
. . . . . . . . . . 150
. . . . . . . . . . 155
irregulares e de
. . . . . . . . . . 158
12 Métodos de iniciação de modelos numéricos
12.1 Iniciação de modelos: método dos mı́nimos quadrados . . . .
12.2 Iniciação de modelos de equações primitivas . . . . . . . . .
12.3 Modelagem hidrodinâmica 2D para uma região geográfica
especı́fica . . . . . . . . . . . . . . . . . . . . . . . . . . . .
164
164
167
171
13 Equações hidrodinâmicas básicas de modelo 3D
173
13.1 Equações básicas de modelo numérico hidrodinâmico 3D . . 173
13.2 Coeficientes de difusão . . . . . . . . . . . . . . . . . . . . . 176
13.3 Modelos de multi-camadas . . . . . . . . . . . . . . . . . . . 177
13.4 Equações básicas de modelo 3D simplificado e inserção de
condições de contorno na superfı́cie . . . . . . . . . . . . . . 179
13.5 Discretização das equações básicas de modelo 3D simplificado 180
14 Modelos numéricos hidrodinâmicos 3D
14.1 Modelo de coordenadas z explı́cito com 2 nı́veis de tempo
14.2 Modelo de coordenadas z explı́cito com 3 nı́veis de tempo
14.3 Discretização de termos com formulação implı́cita . . . .
14.4 Tipos de modelos 3D . . . . . . . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
15 Separação de modos; coordenadas verticais sigma, isopicnais e hı́bridas
15.1 Separação de modos . . . . . . . . . . . . . . . . . . . . . .
15.2 Coordenadas verticais sigma . . . . . . . . . . . . . . . . . .
15.3 Coordenadas verticais isopicnais . . . . . . . . . . . . . . . .
15.4 Coordenadas verticais hı́bridas . . . . . . . . . . . . . . . . .
15.5 Solução espectral na vertical . . . . . . . . . . . . . . . . . .
15.6 Modelagem hidrodinâmica 3D para uma região geográfica
especı́fica . . . . . . . . . . . . . . . . . . . . . . . . . . . .
16 Modelos atuais em Oceanografia
185
185
190
193
194
200
200
202
205
208
209
214
216
iv
Sumário
16.1
16.2
16.3
16.4
Modelos
Modelos
Modelos
Modelos
do transporte de sedimentos . . . . . . . . . . . .
de propagação de ondas de superfı́cie . . . . . . .
de qualidade da água (Eulerianos e Lagrangeanos)
ecológicos . . . . . . . . . . . . . . . . . . . . . .
.
.
.
.
216
217
219
220
17 O método de elementos finitos
225
17.1 Modelos numéricos de elementos finitos 2D (verticalmente
integrados) . . . . . . . . . . . . . . . . . . . . . . . . . . . . 225
17.2 Modelos numéricos de elementos finitos 3D . . . . . . . . . . 234
Referências Bibliográficas
242
v
Lista de Figuras
1.1
1.2
1.3
2.1
2.2
Gráfico da função e derivadas, como calculadas ao final da
execução do programa exp01 01diferencas finitas.m. . . . . . . .
7
Configuração final do resultado do modelo exp01 02adv ordem01.m. 15
Configuração final do resultado do modelo exp01 03adv ordem01.m. 16
Resultado do modelo exp02 01dif explic.m, com sinal retangular
inicial na parte central da grade (em vermelho) e sua configuração
final devido à difusão (em azul). . . . . . . . . . . . . . . . . . .
Resultado do modelo exc02 04adv implic.m, com sinal retangular
inicial na parte central da grade (em vermelho) e sua configuração
final devido à advecção (em azul). . . . . . . . . . . . . . . . . .
3.1
18
24
Solução da equação da advecção de um sinal retangular que se
move numa grade uni-dimensional à velocidade de 2 m/s (posição
inicial do sinal assinalada em vermelho e posição final em azul);
(a) solução analı́tica: ao final de 300 s, o sinal se moveu 600 m; (b)
erro de instabilidade: ao final de apenas 45 s, a solução apresenta
amplitudes errôneas; (c) erro de dispersão: a localização do sinal
advectado (em azul) é diferente da solução analı́tica (em preto),
pois a velocidade de fase da solução é diferente da analı́tica; (d)
existência de modos computacionais (ruı́do) nos resultados do
modelo; e (e) erro de difusão: a forma da solução foi modificada. 27
3.2 Resultado do modelo exp03 01adv exp ord01.m. Advecção de um
sinal retangular inicialmente no centro da grade, com esquema
explı́cito e diferenças finitas de primeira ordem. . . . . . . . . . 38
vi
Lista de Figuras
3.3
3.4
3.5
3.6
3.7
4.1
4.2
4.3
4.4
5.1
5.2
Resultado do modelo exc03 01adv exp ord02.m. Advecção de um
sinal retangular inicialmente no centro da grade, com esquema
explı́cito e diferenças finitas de segunda ordem. . . . . . . . . .
Resultado do modelo exc03 02adv implic ord01.m. Advecção de
um sinal retangular inicialmente no centro da grade, com esquema
implı́cito e diferenças finitas de primeira ordem no tempo e no
espaço. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Resultado do modelo exc03 03adv implic ord02.m. Advecção de
um sinal retangular inicialmente no centro da grade, com esquema
implı́cito e diferenças finitas de segunda ordem no tempo e no
espaço. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Resultado do modelo exc03 04adv semi implic ord01.m. Advecção de um sinal retangular inicialmente no centro da grade,
com esquema semi-implı́cito e diferenças finitas de primeira ordem
no tempo e no espaço. . . . . . . . . . . . . . . . . . . . . . . .
Resultado do modelo exc03 05adv semi implic ord02.m. Advecção de um sinal retangular inicialmente no centro da grade,
com esquema semi-implı́cito e diferenças finitas de segunda ordem
no tempo e no espaço. . . . . . . . . . . . . . . . . . . . . . . .
Resultado do modelo exp04 01adv dif 2d bx ord.m. Advecção e
difusão 2D com esquema de baixa ordem. . . . . . . . . . . . . .
Resultado do modelo exc04 02adv dif 2d alta ord.m. Advecção
e difusão 2D com esquema de alta ordem. . . . . . . . . . . . .
Resultado do modelo exp04 02adv dif dec 2d.m. Simulação do
programa de advecção–difusão–decaimento 2D, aplicado a concentração inicial de poluente. . . . . . . . . . . . . . . . . . . .
Resultado do modelo exc04 03adv dif dec 3d.m. Simulação do
programa de advecção–difusão–decaimento 3D, aplicado a descarga contı́nua de poluente. . . . . . . . . . . . . . . . . . . . .
Resultado do modelo exp05 01adv cc rig.m. Configuração do
sinal senoidal modelado após reflexão no final da grade, devido
ao uso de condição de contorno na forma de parede rı́gida. . . .
Resultado do modelo exc05 01adv cc grad.m. Configuração do
sinal senoidal modelado com o uso de condição de contorno não
gradiente. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
39
40
41
42
43
54
55
56
57
65
66
vii
Lista de Figuras
5.3
5.4
5.5
5.6
5.7
5.8
6.1
6.2
6.3
6.4
6.5
Resultado do modelo exc05 02adv cc extrap.m. Configuração
do sinal senoidal modelado com o uso de condição de contorno
correspondente à extrapolação linear de resultados do modelo
em pontos internos. . . . . . . . . . . . . . . . . . . . . . . . . .
Resultado do modelo exc05 03adv cc lat.m. Configuração do
sinal senoidal modelado com o uso de condição de contorno
correspondente a diferença finita lateral. . . . . . . . . . . . . .
Resultado do modelo exp05 02dif cc rig.m. Configuração do
sinal retangular submetido a difusão com o uso de condição de
contorno na forma de parede rı́gida (nas duas extremidades da
grade). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Resultado do modelo exc05 10dif cc rad.m. Configuração do
sinal retangular submetido a difusão com o uso de condição de
contorno radiacional (nas duas extremidades da grade). . . . . .
Resultado do modelo exc05 10dif cc rad.m. Configuração do
sinal retangular submetido a difusão com o uso de condição de
contorno radiacional (nas duas extremidades da grade). . . . . .
Resultado da simulação do modelo exp05 04disp reg geog.m,
com advecção–difusão–decaimento 2D na Baı́a da Guanabara,
aplicado à descarga contı́nua de uma substância na região. . . .
Resultado do modelo exp06 01filtros1d.m. Função 1D (em azul)
e sua filtragem, com e sem ponderação (em verde e vermelho,
respectivamente). . . . . . . . . . . . . . . . . . . . . . . . . . .
Resultado do modelo exp06 03adv ord01.m. Advecção de um
sinal 1D com solução numérica de 1a ordem, com monitoramento
da conservação (em vermelho, a condição inicial). . . . . . . . .
Resultado do modelo exc06 02adv ord02.m. Advecção de um
sinal 1D com solução numérica de 2a ordem, com monitoramento
da conservação (em vermelho, a condição inicial). . . . . . . . .
Resultado do modelo exc06 03adv dif ord02.m. Advecção de um
sinal 1D com solução numérica de 2a ordem, com monitoramento
da conservação e controle do ruı́do computacional através de
difusão (em vermelho, a condição inicial). . . . . . . . . . . . .
Resultado do modelo exc06 04adv filt continua ord02.m. Advecção de um sinal 1D com solução numérica de 2a ordem, com
monitoramento da conservação e controle do ruı́do computacional
através de filtragem contı́nua (em vermelho, a condição inicial).
67
68
69
70
71
72
77
82
83
84
85
viii
Lista de Figuras
6.6
7.1
7.2
7.3
7.4
Resultado do modelo exc06 05adv filt period ord02.m. Advecção
de um sinal 1D com solução numérica de 2a ordem, com monitoramento da conservação e controle do ruı́do computacional
através de filtragem periódica (em vermelho, a condição inicial).
Resultado do modelo exp07 01adv ordem01.m. Advecção de um
sinal retangular com esquema de 1a ordem no espaço. . . . . . .
Resultado do modelo exc07 01adv ordem03filt5.m. Advecção
de um sinal retangular com esquema de 3a ordem no espaço (e
filtragem periódica dos resultados). . . . . . . . . . . . . . . . .
Resultado do modelo exc07 03adv quick.m. Aplicação do esquema QUICK para a advecção de um sinal 1D retangular. . . .
Resultado do modelo exc07 04adv quick filt3p.m. Aplicação do
esquema QUICK para a advecção de um sinal 1D retangular,
com filtragens periódicas dos resultados do modelo. . . . . . . .
86
89
90
96
97
8.1
Resultado do modelo exp08 01eq helmoltz.m. Solução da equação
de Helmholtz, para função constante no contorno. . . . . . . . . 112
9.1
Esquema de grade alternada para pontos tipo η e u com sua
numeração. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114
Resultado do modelo exp09 01hidr1d 01.m. Distribuições de
elevação e correntes calculadas pelo modelo hidrodinâmico 1D
forçado por oscilação harmônica no ponto inicial da grade, ao
final de 30 minutos de simulação. . . . . . . . . . . . . . . . . . 121
9.2
10.1 Grade 2D alternada, para os pontos tipo η, U e V , definição dos
espaçamentos de grade e numeração dos pontos de grade. . . . . 130
10.2 Resultado do modelo exp10 01hidr2d 01.m. Exemplo de resultado de modelo numérico hidrodinâmico 2D forçado por vento
de Sudoeste, com elevações (esquerda) e correntes (direita) ao
final de 10 minutos de simulação. . . . . . . . . . . . . . . . . . 145
11.1 Grades de Arakawa, dos tipos A (não alternada) e B, C, D e E
(alternadas). Note-se que, em geral, ao utilizar grades alternadas,
se consideram os mesmos ı́ndices (j, k) para cada trio de pontos
η, U e V (como apresentado no Capı́tulo 10). Somente no estudo
comparativo de grades deste Capı́tulo 11 é que serão utilizados ı́ndices
(j, k) para cada posição ao longo dos eixos (x, y). . . . . . . . . . . 151
ix
Lista de Figuras
11.2 Relações de dispersão analı́tica (11.6) e referentes às grades de
Arakawa A–E (equações 11.8 a 11.12). . . . . . . . . . . . . . .
11.3 Aninhamento de grade em sub-região de interesse. . . . . . . . .
11.4 Aninhamento de grades e processamento com influência nos dois
sentidos. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
11.5 Grades principal e aninhada (alternadas) e numeração de seus
pontos, com ı́ndices J, K (maiúsculas) referentes à grade maior e
ı́ndices j, k (minúsculas) referentes à grade menor. . . . . . . . .
11.6 Exemplo de grade inclinada. GP: grade principal; A1: grade
aninhada 1 e A2: grade aninhada 2. . . . . . . . . . . . . . . . .
11.7 Exemplo de grade curvilı́nea. . . . . . . . . . . . . . . . . . . .
11.8 Exemplo de grade irregular (com espaçamento variável ∆xjk ). .
11.9 Grade alternada com contorno terrestre variável no tempo, formada por pontos tipo u (−) e η (+), nesta ordem, com um ı́ndice
j para cada par de pontos u, η. . . . . . . . . . . . . . . . . . .
154
155
156
157
159
160
160
162
12.1 Grade e pontos (×) de coordenadas x0 e y0 com medições f0
a serem usadas no ajuste de mı́nimos quadrados para impor
condições iniciais no processamento de um modelo. . . . . . . . 166
12.2 Resultado do modelo exp12 01hidr2d 01baia.m. Exemplo de
cálculo de modelo numérico hidrodinâmico 2D numa região semifechada, forçado por vento de Sudoeste, com elevações (esquerda)
e correntes (direita) ao final de 10 minutos de simulação. . . . . 170
12.3 Resultado da simulação do modelo exc12 04hidr2d reg geog.m,
com elevações da superfı́cie e correntes médias na vertical, na
Baı́a da Guanabara. . . . . . . . . . . . . . . . . . . . . . . . . 172
13.1 Esquema de modelo de duas camadas. . . . . . . . . . . . . . . 177
13.2 Resultado do modelo exp13 01hidr3d 01.m. Exemplo de cálculo
de modelo numérico hidrodinâmico 3D numa região fechada,
forçado por vento de Sudoeste, com elevações (esquerda) e correntes em vários nı́veis (direita) ao final de 30 minutos de simulação.184
14.1 Parte x, y da grade alternada tri-dimensional x, y, z: + pontos
tipo η, − pontos tipo u e | pontos tipo v (elipses delimitam
pontos com mesmos ı́ndices j, k). . . . . . . . . . . . . . . . . . 196
x
Lista de Figuras
14.2 Parte x, z da grade alternada tri-dimensional x, y, z: + pontos
tipo T S, − pontos tipo u e | pontos tipo w (elipses delimitam
pontos com mesmos ı́ndices j, l). Os coeficientes de difusão
vertical são referentes aos pontos tipo w. . . . . . . . . . . . . . 197
14.3 Parte y, z da grade alternada tri-dimensional x, y, z: + pontos
tipo T S, − pontos tipo v e | pontos tipo w (elipses delimitam
pontos com mesmos ı́ndices k, l). Os coeficientes de difusão
vertical são referentes aos pontos tipo w. . . . . . . . . . . . . . 198
14.4 Grade alternada tri-dimensional x, y, z completa. . . . . . . . . 199
15.1 Modelos de grades com coordenadas verticais sigma, sobre um
monte submarino (a cima, incluindo nı́veis z) e de uma região
costeira para uma área profunda (a baixo). . . . . . . . . . . . . 203
15.2 Esquema de coordenadas hı́bridas, passando de coordenadas
isopicnais (densidades ρ) para coordenadas de nı́vel e para coordenadas sigma. . . . . . . . . . . . . . . . . . . . . . . . . . . . 209
15.3 Resultado da simulação do modelo exc15 01hidr3d reg geog.m,
com elevações da superfı́cie e correntes na superfı́cie, na Baı́a da
Guanabara. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 215
17.1 Grades de modelos de elementos finitos, de baixa resolução (a) e
alta resolução (b). . . . . . . . . . . . . . . . . . . . . . . . . . 226
xi
Lista de Tabelas
1.1
1.2
1.3
4.1
4.2
4.3
4.4
6.1
Evolução dos espaçamentos de grade utilizados em modelos
numéricos em Oceanografia no decorrer das últimas décadas. . .
Diferenças finitas no espaço x e no tempo t, de 1a e 2a ordens,
de 2 e 3 pontos, laterais (avançadas ou retardadas) e centradas.
Nı́veis de tempo e representação da variável f na forma matricial.
Resumo das condições de estabilidade para a equação da advecção.
Resumo das condições de estabilidade para a equação da difusão.
Resumo das condições de estabilidade para a equação do decaimento. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Resumo das condições de estabilidade para a equação do oscilador.
Caracterı́sticas de modelos lineares e não lineares quanto aos
números de onda. . . . . . . . . . . . . . . . . . . . . . . . . . .
3
6
9
45
45
47
48
79
15.1 Esquemas de modelos 3D, 2D e 3D′ . . . . . . . . . . . . . . . . 201
xii
Prefácio
Este livro é produto de muitos anos ministrando cursos de Modelagem
numérica em Oceanografia no Instituto Oceanográfico da Universidade de
São Paulo, desde 1982, ao nı́vel de graduação, pós-graduação e especialização.
Tenho a agradecer imensamente a muitas turmas de alunos que se
dedicaram a este aprendizado e me proporcionaram enorme experiência de
ensino nesta área do conhecimento. Sem as correções, crı́ticas, observações e
sugestões dos alunos esta obra não seria concretizada.
Agradeço também a meus professores e docentes colegas da universidade,
por seus ensinamentos e constantes trocas de informações e experiências.
Agradeço em especial a meus orientadores, Alm. Alberto dos Santos Franco
(in memoriam) e Prof. Dr. Afranio Rubens de Mesquita, e meu colega de
inúmeras pesquisas, Prof. Dr. Ricardo de Camargo.
Agradeço aos funcionários e técnicos do Instituto Oceanográfico da
Universidade de São Paulo, pelo apoio permanente em minhas atividades de
ensino e pesquisa, especialmente a Sra. Vanilde Ferreira de Oliveira.
Um agradecimento especial é dirigido a meus familiares, que tanto me
apoiaram em todas as fases de minha carreira acadêmica, meus pais Cesar
Youssef Harari (‫ )ז ״ל‬e Rachel Harari (‫)ז ״ל‬, minha irmã Camille Harari, meus
sogros Chyia Szajnbok e Dvorah Szajnbok, meus filhos Rachel Szajnbok
Harari e Cesar Szajnbok Harari e, especialmente, minha querida esposa
Gina Szajnbok Harari, autora de diversas cuidadosas revisões do texto.
Finalmente, agradeço a Deus Todo Poderoso, Criador do Universo,
fonte da sabedoria que permite aos homens desenvolverem as ciências e o
conhecimento do Universo.
Joseph Harari
São Paulo, 2015
xiii
Introdução
Este livro tem como objetivo apresentar os conceitos básicos e as técnicas
utilizadas na modelagem numérica em Oceanografia. Pode ser utilizado em
cursos de vários nı́veis—como graduação, pós-graduação ou especialização—
dependendo do grau de aprofundamento determinado pelos professores
responsáveis em cada tópico abordado.
Na parte inicial do livro são apresentados os conceitos básicos em modelagem numérica em Oceanografia, esquemas de diferenças finitas, análise
de erros em modelos numéricos, determinação de condições de estabilidade,
imposição de condições de contorno em modelos e utilização de filtragens de
resultados de modelo.
Ao final da parte introdutória, se encontram esquemas numéricos de alta
precisão, soluções de equações elı́pticas e modelos numéricos hidrodinâmicos
1D (em uma dimensão, ao longo de um canal, por exemplo).
Na parte final do livro, são apresentados numéricos hidrodinâmicos 2D
(em duas dimensões, verticalmente integrados), com a utilização de grades
cobrindo as regiões modeladas e métodos de iniciação de modelos numéricos;
a seguir, são analisados modelos numéricos hidrodinâmicos 3D (em três
dimensões), com vários tipos de solução na vertical.
Como complementação, são apresentados, em linhas gerais, modelos
atuais em Oceanografia, como por exemplo modelos do transporte de sedimentos, modelos de propagação de ondas de superfı́cie, modelos de qualidade
da água, modelos ecológicos e, finalmente, o método de elementos finitos
(em duas e três dimensões).
Em cada capı́tulo, são fornecidos exemplos e exercı́cios com soluções, em
linguagem MATLABr . Os exemplos são identificados como expCAP numeronome.m
e os exercı́cios como excCAP numeronome.m, onde CAP é o número de
xiv
identificação do capı́tulo e numeronome contém o número e a identificação
do exemplo ou exercı́cio.
É altamente recomendável que alunos e usuários do livro tentem resolver
os exercı́cios antes de consultar as soluções apresentadas. Os exemplos
e exercı́cios foram organizados de forma simples e didática, evitando-se
otimizações exageradas e procurando facilitar o seu entendimento pelos
usuários. Também não foram elaborados com a finalidade de uso imediato
em análises oceanográficas, porém podem vir a ser usados com este objetivo,
desde que sejam ajustados os coeficientes e parâmetros de cada programa
para as aplicações desejadas, e que seus resultados sejam devidamente
comparados com medições ou informações independentes.
Finalizando, é solicitado aos leitores deste livro que encaminhem ao autor
correções, aprimoramentos, comentários e sugestões, para que o mesmo possa
ser continuamente aprimorado e possa ser de utilidade a alunos interessados
nas aplicações da modelagem numérica em Oceanografia.
xv
1
Conceitos básicos em Modelagem
Numérica
1.1
Introdução
O método cientı́fico consiste em quatro etapas:
1.
2.
3.
4.
observação de um fenômeno;
medição de parâmetros, de modo a quantificar o fenômeno;
elaboração de teorias para explicar o fenômeno;
reprodução do fenômeno em laboratório ou sua previsão.
Como exemplo da aplicação da metodologia cientı́fica, em meteorologia,
são observados perı́odos de frio em latitudes médias da América do Sul; são
realizadas medições de temperatura, pressão atmosférica, umidade e ventos;
teorias sobre a geração de frentes frias em altas latitudes e sua evolução são
desenvolvidas; e novos perı́odos de frio podem então ser previstos.
Em Oceanografia, a modelagem numérica utiliza medições e teorias sobre
o comportamento do oceano, de modo a possibilitar simulações e previsões
dos processos que nele ocorrem, como a circulação marı́tima, o transporte de
sedimentos, a cadeia alimentar, etc. Em particular, a modelagem numérica
hidrodinâmica utiliza medições de nı́vel do mar, correntes e propriedades
fı́sico–quı́micas da água do mar e resolve numericamente as equações hidrodinâmicas básicas, de modo a reproduzir e prever a circulação marı́tima e a
distribuição de propriedades. De fato, a modelagem da circulação constitui
a base dos demais modelos em Oceanografia, visto que seus resultados são
utilizados na modelagem de ondas, sedimentos, poluentes, etc. . .
Há dois tipos de modelos de circulação marı́tima: os gerais e os especı́ficos.
Os modelos gerais procuram modelar o oceano da maneira mais completa
1
1.1. Introdução
possı́vel, sintetizando conjuntamente diversos processos fı́sicos e adotando
hipóteses simplificadoras e parametrizações dos processos de escala menores.
Um exemplo de um modelo geral é o da circulação geral dos oceanos, que
considera conjuntamente todas as componentes da circulação, devidas a
ventos, marés e variações de densidade. Já os especı́ficos, procuram estudar
fenômenos de forma individual, isolando-os de outros processos (até onde
possı́vel). Um exemplo de modelo especı́fico é um modelo da circulação de
maré, que evidentemente não inclui as componentes de circulação geradas
pelo vento e por variações de densidade.
Uma divisão importante nos modelos em Oceanografia (hidrodinâmicos,
do transporte de sedimentos, da cadeia alimentar etc. . . ) é referente a suas
escalas espaciais:
1. grande escala: para simulações em escala global ou de bacia oceânica;
2. meso-escala: estuda processos em plataformas continentais e
3. pequena escala: cobre regiões costeiras e estuários.
A escala está intimamente relacionada com a resolução de um modelo,
definida pela distancia horizontal entre os pontos de cálculo de uma grade
computacional (δx, “espaçamento de grade”). A evolução da informática,
de certa forma, definiu ao longo dos últimos anos a resolução adotada por
cada escala, como exemplificado na Tabela 1.1.
2
Tabela 1.1: Evolução dos espaçamentos de grade utilizados em modelos
numéricos em Oceanografia no decorrer das últimas décadas.
Espaçamentos de grade [km]
Escalas
Década de 1980
Década de 1990
Década de 2000
Década de 2010
Grande
Meso
Pequena
100
010
001
10,0
01,0
10,1
1,01
0,10
0,01
0,100
0,010
0,001
1.1. Introdução
3
1.2. Método de diferenças finitas
A resolução de um modelo define que fenômenos são efetivamente incluı́dos nas simulações e quais são omitidos ou simplesmente parametrizados,
por possuı́rem distâncias horizontais caracterı́sticas da mesma ordem (ou
menor) que os espaçamentos de grade adotados. É importante notar que
os espaçamentos de grade devem ser muito menores que os comprimentos
de ondas das oscilações modeladas (∆x ≪ L), ou seja, deve haver um
grande número de pontos (e espaçamentos) de grade para representar um
comprimento de onda (no mı́nimo, L = 20 ∆x ou L = 30 ∆x).
Em geral, modelos requerem dados iniciais e de contorno para o seu
processamento, além de dados para calibração e validação; esses dados
podem ser obtidos através de medições in situ ou remotas, ou de resultados
de outros modelos. A calibração consiste na comparação de resultados do
modelo com um conjunto de dados, visando obter os valores dos coeficientes
de fricção e de atenuação mais adequados. Uma vez calibrado um modelo,
a validação consiste na comparação de seus resultados com outro conjunto
de dados, de modo a inferir a qualidade dos resultados das simulações do
modelo.
É evidente que modelos, teorias e observações devem se desenvolver
conjuntamente, afim que haja interação entre os mesmos, com observações
resultando em teorias mais aprimoradas, que conduzem a modelos mais sofisticados, os quais, por sua vez, demandam observações melhores e possibilitam
aprimorar as teorias.
Por fim, vale ressaltar a grande similaridade entre a modelagem em
Meteorologia e Oceanografia, uma vez que as equações que descrevem seus
processos fı́sicos são muito semelhantes (Haltiner, 1971; Haltiner e Williams,
1980; Ramming e Kowalik, 1980; Kowalik e Murty, 1993).
1.2
Método de diferenças finitas
A aplicação do método de diferenças finitas na modelagem numérica
consiste em discretizar o espaço e o tempo, de modo que se possa substituir
derivadas por diferenças finitas (Fortuna, 2000).
A partir da discretização do espaço x com espaçamentos ∆x e ı́ndices
j (x = j ∆x, ∆x > 0), o ponto de partida do método de diferenças finitas
provém da Série de Taylor (1.1) de uma função f (x):
∆x2
∆x3
′′′
f (x + ∆x) = f (x) + f (x)∆x + f (x)
+ f (x)
+ ···
2!
3!
′
′′
(1.1)
4
1.2. Método de diferenças finitas
Isolando a primeira derivada, tem-se:
f ′ (x) =
f (x + ∆x) − f (x)
+R
∆x
(1.2)
∆x
onde o termo de maior ordem no resı́duo R é f ′ (x)
.
2!
Se R for desprezado, a aproximação 1.2 é chamada de diferença finita
avançada, uma vez que ∆x é maior que zero, sendo uma diferença finita de
2 pontos e ordem ∆x, cuja representação é O(∆x). Esta ordem representa o
erro de truncamento das diferenças finitas, o que indica sua precisão.
Contudo, a série de Taylor pode ser expressa também da seguinte forma:
f (x − ∆x) = f (x) − f ′ (x) ∆x + f ′′ (x)
∆x2
∆x3
− f ′′′ (x)
+ ···
2!
3!
(1.3)
que, analogamente ao primeiro caso, conduz a:
f ′ (x) =
f (x) − f (x − ∆x)
+R
∆x
(1.4)
A equação 1.4 representa uma diferença finita retardada, de 2 pontos e
também de ordem de precisão O(∆x).
Ainda é possı́vel subtrair (1.1) de (1.2), resultando em:
f ′ (x) =
f (x + ∆x) − f (x − ∆x)
+R
2∆x
(1.5)
Esta aproximação é uma diferença finita centrada, de 3 pontos e ordem
de precisão O(∆x2 ), sendo portanto superior às outras duas. Embora de
melhor precisão, a diferença centrada não pode ser aplicada no primeiro e no
último ponto da grade, sendo necessárias as diferenças avançada e retardada,
respectivamente.
Também se pode somar (1.1) de (1.2), permitindo determinar a segunda
derivada, em diferença finita centrada, com 3 pontos e ordem de precisão
O(∆x2 ):
f ′′ (x) =
f (x + ∆x) − 2f (x) + f (x − ∆x)
+R
∆x2
(1.6)
Vale lembrar que para o modelo seja de qualidade é necessário que o
espaçamento de grade seja muito menor que o comprimento de onda do
fenômeno a ser estudado.
5
1.2. Método de diferenças finitas
Tabela 1.2: Diferenças finitas no espaço x e no tempo t, de 1a e 2a ordens,
de 2 e 3 pontos, laterais (avançadas ou retardadas) e centradas.
Espaço
Tempo
∂f
fj+1 − fj
=
∂x
∆x
∂f
f n+1 − f n
=
∂t
∆t
∂f
fj − fj−1
=
∂x
∆x
∂f
f n − f n−1
=
∂t
∆t
∂f
fj+1 − fj−1
=
∂x
2∆x
∂f
f n+1 − f n−1
=
∂t
2∆t
fj+1 − 2fj fj−1
∂ 2f
=
2
∂x
∆x2
f n+1 − 2f n + f n−1
∂f
=
∂t2
∆t2
1a ordem
2 pontos
avançada
1a ordem
2 pontos
retardada
2a ordem
3 pontos
centrada
2a ordem
3 pontos
centrada
Expressões similares a (1.2), (1.4), (1.5) e (1.6), referentes à coordenada
x, podem ser estabelecidas para derivadas em relação às coordenadas y e z
e ao tempo t. Para espaçamentos ∆x, ∆y, ∆z, ∆t e ı́ndices j, k, l, n, são
consideradas discretizações da forma x = j ∆x , y = k ∆y , z = l ∆z e
t = n ∆t, respectivamente. Uma função de x, y, z, t será representada por:
n
f = f (x, y, z, t) = f (j∆x, k∆y, l∆z, n∆t) = fj,k,l
(1.7)
de modo que as diferenças finitas, escritas em termos de ı́ndices, para x e t,
se encontram na Tabela 1.2.
Exemplo 1.1: exp01 01diferencas finitas.m
programa exp01 01diferencas finitas.m, para cálculo de primeiras e
segundas derivadas da função f (x) = sen(x), através de fórmulas
analı́ticas e de diferenças finitas de 1a e 2a ordem, utilizando as
relações (1.2), (1.4), (1.5) e (1.6), e incluindo estatı́sticas dos erros
(ver Figura 1.1).
6
1.3. Diferenças finitas de alta ordem
Função (azul) e diferenças finitas
avançada (vermelho), retardada (verde), centrada (preto)
f,
df
alta ordem
dx
1
0,5
0
−0,5
−1
0
0,5
1
1,5
x
2
2,5
3
Figura 1.1: Gráfico da função e derivadas, como calculadas ao final da
execução do programa exp01 01diferencas finitas.m.
1.3
Diferenças finitas de alta ordem
Diferenças finitas laterais de segunda ordem podem ser obtidas através
da utilização das expressões de Séries de Taylor (1.1) e (1.3), levando em
conta aproximações da segunda derivada (1.6) centradas nos pontos (x + ∆x)
e (x − ∆x); resultam então, para as diferenças finitas avançada e retardada
de segunda ordem (Mitchell e Griffiths, 1980),
∂f
−3fj + 4fj+1 − fj+2
=
+R
∂x
2∆x
3fj − 4fj−1 + fj−2
∂f
=
+R
∂x
2∆x
(1.8)
(1.9)
Tomando a média dessas duas expressões resulta uma diferença finita
centrada de terceira ordem:
fj+1 − fj−1 fj+2 − fj−2
∂f
=2
−
+R
∂x
2∆x
4∆x
(1.10)
7
1.4. Equação da advecção
Exercı́cio 1.1: exc01 01diferencas finitas alta ordem.m
desenvolver o programa exc01 01diferencas finitas alta ordem.m, para
calcular as primeiras derivadas da função f (x) = sen(x), através de
fórmulas analı́ticas e de diferenças finitas de 2a e 3a ordem, utilizando
as relações (1.8), (1.9) e (1.10), e incluindo estatı́sticas dos erros.
Exercı́cio 1.2: exc01 02diferencas finitas qq f.m
desenvolver o programa exc01 02diferencas finitas qq f.m, para calcular as primeiras derivadas da função f (x) = x sen(x) exp(−x2 ), através
de fórmulas de diferenças finitas de 1a , 2a e 3a ordem, utilizando as
relações (1.6), (1.5), (1.6), (1.8), (1.9) e (1.10).
1.4
Equação da advecção
Uma primeira abordagem do método das diferenças finitas pode ser feita
em uma equação linear simples com coeficientes constantes, comumente
chamada de equação da advecção, que representa o transporte de uma
propriedade:
∂f
∂f
+c
=0
(1.11)
∂t
∂x
onde t é o tempo, x o espaço, f a propriedade advectada e c, diferente de
zero, é a velocidade, ou ainda, a velocidade de fase, no caso de se considerar
ondas. A equação (1.11) é similar à versão linearizada dos primeiros termos
da equação do movimento, da equação da termodinâmica e da equação
da vorticidade, entre outras. A equação da advecção considerada pode
ser aproximada pelas diferenças descritas acima, sendo que, inicialmente,
será utilizada a diferença finita avançada no tempo e retardada no espaço
(portanto, de 1a ordem):
n
fjn+1 − fjn
fjn − fj−1
+c
=0
∆t
∆x
(1.12)
cuja fórmula de recorrência é:
fjn+1 = fjn − c
∆t n
n
fj − fj−1
∆x
(1.13)
8
1.5. Estrutura de um modelo
Tabela 1.3: Nı́veis de tempo e representação da variável f na forma matricial.
Nı́veis de tempo
Variáveis
Matrizes
Instante
(n − 1)
(n)
(n + 1)
f n−1
fn
f n+1
fant
fatu
fren
“anterior”
“atual”
“renovado”
onde ∆x corresponde ao espaçamento de grade e ∆t, ao passo no tempo. A
seguir, a Tabela 1.3 apresenta a notação matricial (dependente do espaço)
da variável f em função dos nı́veis de tempo n − 1, n, n + 1. Ao programar
a fórmula de recorrência (1.13), a variável f é considerada em dois nı́veis
de tempo, sendo portanto representada no espaço por duas matrizes, fren e
fatu (ver Tabela 1.3).
Portanto, a equação de recorrência para a equação da advecção (1.13),
na forma matricial, pode ser escrita como:
f ren (j) = f atu (j) − c
∆t
(f atu (j) − f atu (j − 1))
∆x
(1.14)
A solução da equação da advecção (1.11) por diferenças finitas centradas
no tempo e no espaço, de 2a ordem portanto, conduz a:
n
n
fjn+1 − fjn−1
fj+1
− fj−1
+c
=0
2∆t
2∆x
(1.15)
cuja relação de recorrência é:
∆t n
n
fj+1 − fj−1
∆x
(1.16)
∆t
(f atu(j + 1) − f atu(j − 1))
∆x
(1.17)
fjn+1 = fjn−1 − c
e sua forma matricial corresponde a:
f ren(j) = f ant(j) − c
1.5
Estrutura de um modelo
Primeiramente, vale lembrar que as variáveis em um modelo são usualmente representadas como matrizes espaciais. Além disso, uma variável
9
1.6. Método de elementos finitos
em instantes distintos de tempo é representada através de uma matriz para
cada instante, a qual é renovada a cada passo de tempo. Desta forma, nas
equações discretizadas segundo a forma matricial, se considera que: o termo
fjn−1 corresponde à matriz f ant (variável no instante de tempo anterior), o
termo fjn corresponde a f atu (variável no instante atual) e fjn+1 corresponde
à matriz f ren (variável no instante renovado).
Assim, um modelo pode ser estruturado segundo as seguintes etapas:
1.
2.
3.
4.
5.
6.
7.
8.
9.
Definição dos parâmetros do modelo (∆x, ∆t, c, etc.)
Estabelecimento da grade computacional
Cálculo de constantes iniciais
Imposição de condições iniciais, por exemplo: fjn=1 e fjn=2
Inı́cio do laço (loop) no tempo
Especificação de condições de contorno, quando houver
Evolução da solução no tempo com a fórmula de recorrência
Gravação de resultados (gráficos, arquivos, etc.)
Transferência de variáveis no tempo (f ant = f atu e f atu = f ren,
nesta ordem)
10. Retorno ao inı́cio do laço no tempo, até atingir o número total de
iterações
11. Final do processamento
1.6
Método de elementos finitos
Quanto à solução numérica das equações que compõem um modelo, há
duas formulações: 1) no método de diferenças finitas, acima apresentado, as
derivadas são substituı́das por diferenças finitas, para pontos de grade no
espaço e nı́veis de tempo (Ferziger e Peric, 1997; Strickwerda, 1989); 2) no
segundo método, o de elementos finitos, uma solução é expressa em termos
de uma expansão com funções base que levam em conta sua dependência
espacial (pré-estabelecidas), tendo cada função base um coeficiente que
depende do tempo (a determinar). Ao substituir as expansões nas equações
do modelo, resultam equações diferenciais ordinárias para os coeficientes;
esses coeficientes, por sua vez, evoluem no tempo através de cálculos por
diferenças finitas (Zienkiewicz e Morgan, 1984; Reddy, 1984).
10
1.6. Método de elementos finitos
Um exemplo da aplicação do método de elementos finitos é dado para a
equação da advecção linear 1D aplicada num canal de comprimento L:
∂f
∂f
+c
=0
∂t
∂x
(1.18)
É considerada uma expansão da solução f (x, t) com funções base dependentes do espaço pré-escolhidas fr (x) e coeficientes dependentes do tempo a
determinar φr (t); sendo r o ı́ndice de cada termo da expansão:
f (x, t) =
∞
X
φr (t) fr (x)
(1.19)
r=1
Ao substituir (1.19) em (1.18) resulta:
∂
∞
X
∂
φr (t) fr (x)
r=1
∂t
+c
∞
X
φr (t) fr (x)
r=1
∂x
=0
(1.20)
Levando em conta a dependência das variáveis no tempo e no espaço
resulta:
∞
∞
X
X
∂fr
∂φr
+c
=0
(1.21)
φr
fr
∂t
∂x
r=1
r=1
Multiplicando esta equação por fs e integrando no domı́nio [0, L]:
L
∞ Z
X
r=1 0
L
∞ Z
X
∂fr
∂φr
fs fr dx
+c
dx φr = 0
fs
∂t
∂x
r=1
(1.22)
0
que pode ser escrito na forma matricial (para um número limitado de termos
da expansão, suficientemente grande)

 L
 L

Z
Z
 fs fr dx ∂φr + c  fs ∂fr dx [φr ] = 0
(1.23)
∂t
∂x
0
0
−1  L

 L
Z
Z
∂fr 
∂φr
(1.24)
dx [φr ] = 0
+ c  fs fr dx  fs
∂t
∂x
0
0
11
1.6. Método de elementos finitos
a qual é resolvida por diferenças finitas:
φn+1
r
φn−1
r

ZL
−1 
ZL

∂fr  n
dx [φr ] = 0
∂x
0
0

−1  L
 L
Z
Z
n+1 n−1 ∂fr  n
dx [φr ] = 0
φr
− φr
+ 2∆t c  fs fr dx  fs
∂x
0
0
−1  L
 L

Z
Z
n+1 n−1 ∂fr  n
− 2∆t c  fs fr dx  fs
= φr
φr
dx [φr ] = 0
∂x
−
2∆t
+ c
0
fs fr dx

fs
(1.25)
(1.26)
(1.27)
0
Como as funções base fr são pré-escolhidas, as integrais que contém essas
funções e suas derivadas em (1.27) podem ser calculadas antes da integração
desta equação no tempo. Uma vez renovados os valores dos coeficientes φr ,
a expansão (1.19) é utilizada para renovar a solução (no espaço e no tempo).
Vários tipos de funções podem ser considerados para as funções base, como
por exemplo funções trigonométricas da forma:
fr (x) = cos [(r − 1) πx/L]
(1.28)
Finalizando esta seção, é importante notar que a vantagem do método
de diferenças finitas está na formulação matemática simplificada, seja nos
procedimentos de avanço no espaço e no tempo dos esquemas explı́citos,
seja na transformação de equações de 3 pontos em equações de 2 pontos
dos esquemas implı́citos e semi-implı́citos. Por outro lado, a vantagem
do método de elementos finitos se encontra na flexibilidade da escolha de
pontos de grade no domı́nio espacial do modelo, a qual pode contemplar alta
densidade de pontos em sub-áreas de interesse. Por outro lado, esquemas de
diferenças finitas podem também adotar grades com espaçamentos variáveis.
Na presente publicação, serão apresentadas as técnicas de modelagem
baseadas no método de diferenças finitas; a utilização de elementos finitos
será restrita aos capı́tulos finais, em modelos numéricos hidrodinâmicos
bi-dimensionais (para as duas coordenadas espaciais) e tri-dimensionais
(somente para a coordenada vertical ou para as três coordenadas espaciais
conjuntamente).
12
1.7. Validação dos resultados de um modelo
1.7
Validação dos resultados de um modelo
Um dos aspectos mais importantes da modelagem se encontra na sua
validação, ou seja, na comparação de seus resultados com outro conjunto
de dados, de modo a inferir a qualidade dos resultados das simulações do
modelo (Cecı́lio, 2006). Há vários métodos estatı́sticos para esta finalidade,
que serão abordados a seguir, considerando X n resultados do modelo e Y n
observações em N instantes de tempo, que possuem valores médios Xmed e
Y med, respectivamente (n é o ı́ndice de tempo, n = 1 : N ).
O primeiro parâmetro comparativo é o coeficiente de correlação linear:
X
(X n − Xmed) (Y n − Y med)
corr = X q n
Xq
2
n
(X − Xmed) ·
(Y n − Y med)2
n
(1.29)
n
Este coeficiente varia entre −1 e 1, com os seguintes significados:
corr = 1: há uma correlação perfeita positiva entre as séries Xn e Yn;
corr = −1: há uma correlação negativa perfeita entre as duas séries, isto
é, se uma aumenta, a outra diminui (e vice versa); e
corr = 0: as duas séries não dependem linearmente uma da outra; no
entanto, pode existir uma dependência não linear entre elas e, portanto, a
comparação deve ser investigada por outros meios.
O segundo parâmetro comparativo é o “erro absoluto médio” (EAM),
calculado como:
X
|X n − Y n |
(1.30)
EAM =
N
O terceiro parâmetro é o “erro absoluto médio relativo à média” (EAMR),
dado por:
EAM
EAM R =
(1.31)
Y med
Outra forma muito conveniente de verificar a concordância entre resultados de um modelo X n com observações Y n (ou outras informações
independentes) foi proposta por Wilmott (1981), que calcula um parâmetro
“Skill”, o qual é mais próximo a 1 quanto melhor o ajuste entre as séries e
mais próximo a zero quanto maior o desajuste, sendo o cálculo baseado na
13
1.7. Validação dos resultados de um modelo
seguinte expressão:
S =1− X
X
(|X n − Y n |)2
(|X n − Y med| + |Y n − Y med|)2
(1.32)
Exemplo 1.2: exp01 02adv ordem01.m
o programa exp01 02adv ordem01.m representa a advecção de um
sinal com variação senoidal no tempo no ponto inicial da grade, com
diferenças finitas de primeira ordem (1.12).
Exercı́cio 1.3: exc01 03adv ordem02.m
desenvolver o programa exc01 03adv ordem02.m, análogo ao do exemplo exp01 02adv ordem01.m, porém com solução de 2a ordem, com
diferenças finitas centradas no tempo e no espaço (1.15).
Comparando os resultados dos dois programas de advecção (Figuras 1.2
e 1.3), se nota que, para o mesmo passo de tempo, o de 1a ordem atenua a
amplitude da oscilação no interior da grade, o que não ocorre com a solução
de 2a ordem.
14
1.7. Validação dos resultados de um modelo
Advecção de sinal senoidal na borda (1a ordem)
tempo 360 segundos
0,6
0,4
[m]
0,2
0
−0,2
−0,4
−0,6
0
200
400
600 800 1000 1200 1400 1600 1800
Distância na grade [m]
final
do
Figura 1.2: Configuração
exp01 02adv ordem01.m.
resultado
do
modelo
15
1.7. Validação dos resultados de um modelo
Advecção de sinal senoidal na borda (2a ordem)
tempo 360 segundos
0,6
0,4
[m]
0,2
0
−0,2
−0,4
−0,6
0
200
400
600 800 1000 1200 1400 1600 1800
Distância na grade [m]
final
do
Figura 1.3: Configuração
exp01 03adv ordem01.m.
resultado
do
modelo
16