High-order CE/SE scheme for dam-break flow simulation

<!--A. Rodrigues de Melo<sup>1,2</sup>, L.M. Gramani<sup>1</sup>, E. Kaviski<sup>1</sup>

<sup>1</sup> <small>Programa de Pós-Graduação em Métodos Numéricos em Engenharia, UFPR,</small><br />
<small>Centro Politécnico, Curitiba, Paraná, Brasil</small><br />
<sup>2</sup> <small>Instituto Federal de Educação</small><br />
<small>Ciência e Tecnologia Catarinense, Araquari, Santa Catarina, Brasil
</small>
-->
==Resumo==

Este artigo apresenta um novo esquema explícito para a solução das equações de águas rasas em uma e duas  dimensões, desenvolvido a partir do método dos elementos de conservação e elementos de solução espaço-tempo,  aqui abreviado por método CE/SE. As funções de base utilizadas são expansões de Taylor de segunda  ordem no tempo e no espaço. Esse aumento na ordem das funções de aproximação produz o aumento  no número de variáveis de marcha no tempo, por isso, além das variáveis fluxo e suas inclinações, também são  incógnitas no presente esquema suas derivadas espaciais de segunda ordem. Um processo iterativo para o  cálculo das derivadas de primeira e segunda ordem é formulado para problemas com choques e descontinuidades.  Experimentos computacionais demonstram acurácia de terceira ordem. Os problemas de ruptura de  barragem unidimensional e bidimensional considerados validam a acurácia e robustez do esquema.

''Palavras-Chave:'' Leis de conservação,  Equações de águas rasas,  Volume de controle espaço-tempo,  Alta ordem de acurácia.

==Abstract==

This paper presents a new explicit scheme for the solution of shallow water equations in one and two  space dimensions, developed from the space-time conservation element and solution element (CE/SE) method. The  basis functions used are second-order Taylor expansions in time and space. This increase in the order of the  approximation functions produces an increase in the number of unknowns in the scheme, therefore, besides the  flow variables and their slopes, their second-order partial derivatives are also unknown in the present  scheme. An iterative process for the calculation of the first and second order derivatives is formulated for  problems with shocks and discontinuities. Computational experiments demonstrate third-order accuracy. The  one-dimensional and two-dimensional dam-break problems presented validate the accuracy and robustness of this  scheme.



==1 Introdução==

As equações de águas rasas constituem um dos modelos mais comumente usados na análise de fluxo de água em  rios ou áreas costais <span id='citeF-1'></span>[[#cite-1|[1]]], sobretudo no que tange à simulações de fluxo de ruptura de barragem  <span id='citeF-2'></span><span id='citeF-3'></span><span id='citeF-4'></span><span id='citeF-5'></span>[[#cite-2|[2,3,4,5]]], parte essencial do projeto e avaliação da  segurança de barragens, controle de inundação de rios, mitigação de desastres de bacias hidrográficas, etc   <span id='citeF-6'></span>[[#cite-6|[6]]]. Além da importância prática em hidráulica e engenharia costal, fornecem um ótimo modelo  matemático para equações diferenciais hiperbólicas não lineares que podem ter soluções como ondas de choque  <span id='citeF-7'></span>[[#cite-7|[7]]] e, por isso, também são frequentemente utilizadas como problema de referência para novos  esquemas de previsão numérica <span id='citeF-8'></span><span id='citeF-9'></span>[[#cite-8|[8,9]]].

Extensivos estudos numéricos foram feitos no sentido de analisar os fenômenos governados por essas equações e  diversos métodos foram desenvolvidos, como de volumes finitos <span id='citeF-10'></span><span id='citeF-11'></span><span id='citeF-12'></span>[[#cite-10|[10,11,12]]], elementos finitos <span id='citeF-13'></span><span id='citeF-14'></span><span id='citeF-15'></span><span id='citeF-16'></span>[[#cite-13|[13,14,15,16]]],  diferenças finitas <span id='citeF-17'></span><span id='citeF-18'></span>[[#cite-17|[17,18]]] e outros <span id='citeF-19'></span><span id='citeF-20'></span>[[#cite-19|[19,20]]].

Recentemente, um novo método baseado em volumes de controle, denominado método dos elementos de  conservação e elementos de solução espaço-tempo, do inglês, ''space-time conservation element  and solution element method (CE/SE)'', foi proposto por Chang e To (1991)  <span id='citeF-21'></span>[[#cite-21|[21]]] e aplicações do mesmo em problemas de ruptura de barragem logo foram feitas  <span id='citeF-22'></span><span id='citeF-23'></span><span id='citeF-24'></span>[[#cite-22|[22,23,24]]]. Este método foi desenvolvido com o objetivo de obter soluções  numéricas de leis de conservação e possui características importantes <span id='citeF-25'></span><span id='citeF-26'></span><span id='citeF-27'></span><span id='citeF-28'></span>[[#cite-25|[25,26,27,28]]]: sua construção combina informações de ambas as formas diferencial e integral das leis  de conservação; o conjunto de variáveis de marcha é formado pelas variáveis fluxo e suas derivadas espaciais;  não são utilizadas técnicas de interpolação ou extrapolação nos valores da malha; não são utilizadas técnicas  baseadas nas características, o que torna o método de fácil extensão à várias dimensões, bem como aplicável à  equações não hiperbólicas, como as equações de Navier-Stokes, por exemplo; possui uma simples notação  estêncil, o que facilita a programação. É caracterizado, também, por sua versatilidade, tendo em vista que  pode ser adaptado à malhas não uniformes e não estruturadas <span id='citeF-29'></span><span id='citeF-30'></span>[[#cite-29|[29,30]]], bem como pode ser  estendido à altas ordens <span id='citeF-31'></span><span id='citeF-32'></span>[[#cite-31|[31,32]]], razão pela qual possui aplicações nas mais diversas  áreas <span id='citeF-33'></span><span id='citeF-34'></span><span id='citeF-35'></span><span id='citeF-36'></span>[[#cite-33|[33,34,35,36]]].

O objetivo deste trabalho é apresentar a construção de um novo esquema CE/SE explícito de alta ordem para as  equações de águas rasas e verificar sua eficiência em problemas de ruptura de barragem. Para isso, este  artigo está organizado da seguinte maneira: na seção [[#2 Equações Governantes e Método CE/SE|2]] é apresentado o modelo matemático a ser  estudado (subseção [[#2.1 Equações de Águas Rasas|2.1]]), bem como é desenvolvido o esquema numérico discreto  (subseção [[#2.2 Esquema CE/SE|2.2]]); na seção [[#3 Exemplos Numéricos|3]], experimentos numéricos unidimensionais e bidimensionais  são considerados com a finalidade de validar o método; as conclusões são tecidas junto à seção [[#4 Conclusões|4]] e,  ao final, na seção [[#5 Referências|5]], dispõe-se a lista de referências.

==2 Equações Governantes e Método CE/SE==

===2.1 Equações de Águas Rasas===

As equações de águas rasas, também conhecidas como equações de Saint Venant, possuem a seguinte forma  vetorial conservativa

<span id="eq-1"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\displaystyle \frac{\partial \textbf{Q}}{\partial t}+\frac{\partial  \textbf{F(Q)}}{\partial x}+\frac{\partial \textbf{G(Q)}}{\partial y} =  \textbf{S(Q)},  </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (1)
|}

em que

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\textbf{Q} = \left[\begin{array}{c}q_1 \\ q_2 \\ q_3 \end{array}\right] = \left[\begin{array}{c}h \\ hu \\ hv \end{array}\right], \quad  \textbf{F} = \left[\begin{array}{c}f_1 \\ f_2 \\ f_3 \end{array}\right] =  \left[ \begin{array}{c}hu \\ hu^2+gh^2/2 \\ huv  \end{array} \right] =  \left[ \begin{array}{c}q_2 \\ q_2^2/q_1+gq_1^2/2 \\ q_2q_3/q_1  \end{array} \right], </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (2)
|}

<span id="eq-3"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\textbf{G} = \left[\begin{array}{c}g_1 \\ g_2 \\ g_3 \end{array}\right] =  \left[ \begin{array}{c}hv \\ huv \\ hv^2+gh^2/2  \end{array} \right] =  \left[ \begin{array}{c}q_3 \\  q_2q_3/q_1\\ q_3^2/q_1+gq_1^2/2  \end{array} \right]\quad \textrm{e} \quad  \textbf{S} = \left[\begin{array}{c}S_1 \\ S_2 \\ S_3 \end{array}\right] =  \left[ \begin{array}{c}0 \\ gq_1(S_{0x}-S_{fx}) \\ gq_1(S_{0y}-S_{fy})  \end{array} \right],  </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (3)
|}

sendo <math display="inline">h</math> a profundidade do fluxo (m); <math display="inline">u</math> e <math display="inline">v</math> são as velocidades médias  do escoamento (m/s) nas direções <math display="inline">x</math> e <math display="inline">y</math>, respectivamente; <math display="inline">g</math> é a aceleração da gravidade (<math display="inline">9.812</math>m/s<math display="inline">^2</math>); <math display="inline">S_{0{x}} = \partial Z_b/\partial x</math>  e <math display="inline">S_{0{y}} = \partial Z_b/ \partial y</math> são as inclinações do fundo do  canal cuja topografia é dada por <math display="inline">Z_b</math>; <math display="inline">S_{f{x}}</math>  e <math display="inline">S_{f{y}}</math> são os termos de resistência ao escoamento,  nas direções <math display="inline">x</math> e <math display="inline">y</math>, respectivamente, definidos por

<span id="eq-4"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>S_{fx} = \displaystyle \frac{n^2u\sqrt{u^2+v^2}}{h^{4/3}} =  \frac{n^2q_2\sqrt{q_2^2+q_3^2}}{q_1^{10/3}}  \quad \textrm{e} \quad S_{fy} = \displaystyle \frac{n^2v\sqrt{u^2+v^2}}{h^{4/3}}  = \frac{n^2q_3\sqrt{q_2^2+q_3^2}}{q_1^{10/3}}, </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (4)
|}

em que <math display="inline">n</math> (s/m<math display="inline">^{1/3}</math>) é o coeficiente de fricção de Manning.

===2.2 Esquema CE/SE===

Considere, por simplicidade, a Eq. ([[#eq-1|1]]) como

<span id="eq-5"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\displaystyle \frac{\partial q_m}{\partial{t}}+\frac{\partial f_m}{\partial x}+\frac{\partial g_m}{\partial  y} = S_m, \quad m = 1,2,3 </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (5)
|}

então, pelo teorema da divergência no espaço <math display="inline">\mathbb{R}^3</math>, tem-se que a Eq. ([[#eq-5|5]])  representa a forma diferencial da lei integral de conservação

<span id="eq-6"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\oint _{S(V)}\textbf{H}_m\cdot d\textbf{s} = \int _V S_mdV, </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (6)
|}

em que <math display="inline">\textbf{H}_m = (f_m, g_m,q_m)</math>, <math display="inline">m = 1,2,3</math> e <math display="inline">S(V)</math> representa o contorno de uma região  espaço-tempo <math display="inline">V \subset \mathbb{R}^3</math>.

No método CE/SE, um elemento de conservação (CE) é uma região espaço-tempo em que a conservação do fluxo, Eq.  ([[#eq-6|6]]), é forçada, enquanto que um elemento de solução (SE) é uma região espaço-tempo, normalmente  distinta, em que as variáveis fluxo são supostamente suaves e a Eq. ([[#eq-5|5]]) é válida. Para definir estes  dois importante objetos, considere primeiramente uma malha no plano <math display="inline">x</math>-<math display="inline">y</math>, conforme Fig. [[#img-1a|1a]]. Existem  dois grupos de pontos, marcados por  círculos e cruzes, que representam nós da malha em dois níveis de tempo  diferentes. A cada ponto <math display="inline">(i,j,k)</math> da malha associa-se um CE e um SE. O CE é definido como o quadrilátero <math display="inline">EFGHE'F'G'H'</math> e o SE é a união do quadrilátero <math display="inline">P''Q''R''S''P'Q'R'S'</math> e o polígono <math display="inline">EFGH</math> (veja  Fig. [[#img-1b|1b]]).

<div id='img-1a'></div>
<div id='img-1b'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|style="text-align: center;padding-top:10px;"|[[Image:draft_Aparicio Nogué_501194349-CE_2D_2.png|240px|]]
|style="text-align: center;padding-top:10px;"|[[Image:draft_Aparicio Nogué_501194349-CE_2D.png|300px|a) Pontos da malha representativa no plano x-y e (b) as definições dos  elementos de solução SE(i,j,k) e elementos de conservação CE(i,j,k) no ponto (i,j,k) da malha  [2012_Zhang, 1999_Zhang].]]
|-
|style="text-align: center;font-size:75%;"|(a)
|style="text-align: center;font-size:75%;"|(b)
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="text-align: center;padding-top:10px;" | '''Figura 1:''' (a) Pontos da malha representativa no plano <math>x</math>-<math>y</math> e (b) as definições dos  elementos de solução SE<math>(i,j,k)</math> e elementos de conservação CE<math>(i,j,k)</math> no ponto <math>(i,j,k)</math> da malha  <span id='citeF-23'></span>[[#cite-23|[23,37]]].
|}

Para todo <math display="inline">(x,y,t) \in \mbox{SE}(i,j,k)</math>, aproxima-se <math display="inline">q_m</math> por um polinômio de Taylor de segunda ordem

<span id="eq-7"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>q_m^*(x,y,t;i,j,k) =(q_m)_{i,j}^k+[(q_m)_x]_{i,j}^k(x-x_i)+</math><math>[(q_m)_y]_{i,j}^k(y-y_j)+[(q_m)_t]_{i,j}^k(t-t_k)+[(q_m)_{xt}]_{i,j}  ^k(x-x_i)(t-t_k) </math>
|-
| style="text-align: center;" | <math>   +[(q_m)_{yt}]_{i,j}^k(y-y_j)(t-t_k)+[(q_m)_{xy}]_{i,j}^k(x-x_i)(y-y_j)+\frac{1}{2}[(q_m)_{xx}]_{i,j} ^k(x-x_i)^2 </math>
|-
| style="text-align: center;" | <math>    +\frac{1}{2}[(q_m)_{yy}]_{i,j}^k(y-y_j)^2+\frac{1}{2}[(q_m)_{tt}]_{i,j}^k(t-t_k)^2, \quad m = 1,2,3. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (7)
|}

Aproximações análogas são feitas sobre as funções <math display="inline">f_m</math> e <math display="inline">g_m</math>, <math display="inline">m = 1,2,3</math>. O termo fonte, por outro lado,  é aproximado por uma ordem a menos:

<span id="eq-8"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>S_m^*(x,y,t;i,j,k) =  (S_m)_{i,j}^k+</math><math>[(S_m)_x]_{i,j}^k(x-x_i)+[(S_m)_y]_{i,j}^k(y-y_j)+[(S_m)_t]_{i,j}^k(t-t_k), \quad m = 1,2,3. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (8)
|}

Dessa forma, a Eq. ([[#eq-5|5]]) pode ser  aproximada no SE<math display="inline">(i,j,k)</math> por div <math display="inline">\textbf{H}_m^* = S_m^*</math>, em  que <math display="inline">\textbf{H}_m^* = \left(f_m^*, g_m^*, q_m^*\right)</math>, <math display="inline">m = 1,2,3</math>, ao mesmo tempo em que a Eq.  ([[#eq-6|6]]) é aproximada no CE<math display="inline">(i,j,k)</math> por

<span id="eq-9"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\oint _{S(\textrm{CE}(i,j,k))}\textbf{H}_m^*\cdot d\textbf{s} =  \int _{\textrm{CE}(i,j,k)}  S_m^* dV,  \quad  m = 1,2,3. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (9)
|}

====2.2.1 Avaliação de (q<sub>m</sub>)<sub>i,j</sub><sup>k</sup>====

Substituindo as funções <math display="inline">q_m^*, f_m^*, g_m^*</math> e <math display="inline">S_m^*</math> em ([[#eq-9|9]]), após diversas simplificações,  obtém-se o esquema de avanço no tempo para a variável <math display="inline">(q_m)_{i,j}^k</math>:

<span id="eq-10"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>(q_m)_{i,j}^k + \displaystyle \frac{\Delta x^2}{48}[(q_m)_{xx}]_{i,j}^k+\frac{\Delta  y^2}{48}[(q_m)_{yy}]_{i,j}^k = [P_m^{y}(\Delta x, f_m,g_m)]_{i-1/2,j}^{k-1/2}+[P_m^{y}(-\Delta  x, f_m, g_m)]_{i+1/2,j}^{k-1/2}+ </math>
|-
| style="text-align: center;" | <math>   [P_m^{x}(\Delta y, g_m, f_m)]_{i,j-1/2}^{k-1/2}+[P_m^{x}(-\Delta y, g_m, f_m)]_{i,j+1/2}^{k-1/2}, \quad m  = 1,2,3 \quad \quad  </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (10)
|}

em que

<span id="eq-11"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>P_m^{\zeta }(\Delta \xi ,\phi ,\psi ) = \displaystyle  \frac{1}{4}\left\{ q_m+\frac{\Delta  \xi }{4}(q_m)_{\xi }+\frac{\Delta t}{4}\left[\phi _{\xi }-\psi _{\zeta }\right]+\frac{2\Delta t}{\Delta  \xi }\left[\phi{+\displaystyle}\frac{\Delta t}{4}\phi _t \right]+\right.</math><math> \frac{\Delta t^2}{16}\left[\phi _{\xi  t}-\psi _{\zeta t} \right]+ \frac{\Delta t}{48\Delta \xi }\left[\Delta \zeta ^2\phi _{\zeta  \zeta }+\Delta \xi ^2\phi _{\xi \xi }\right]</math>
|-
| style="text-align: center;" | <math>   +\frac{\Delta t}{12\Delta \xi }\left[\Delta t^2\phi _{tt}-\frac{\Delta  \xi ^2}{2}\psi _{\zeta \xi }\right]+ \frac{\Delta \zeta ^2}{192}(q_m)_{\zeta \zeta }+\frac{7\Delta  \xi ^2}{192}(q_m)_{\xi \xi }+</math><math>\left.\frac{\Delta t}{2}\left[S_m+\frac{\Delta t}{4}(S_m)_t+\frac{\Delta  \xi }{4}(S_m)_{\xi }\right]\right\}.  </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (11)
|}

É importante notar que as Equações ([[#eq-10|10]]) e ([[#eq-11|11]]) dependem apenas das incógnitas <math display="inline">q_m,  (q_m)_x, (q_m)_y, (q_m)_{xx}, (q_m)_{yy}</math> e <math display="inline">(q_m)_{xy}</math> do tempo <math display="inline">t_{k-1/2}</math>, para <math display="inline">m = 1,2,3</math>  pois, pelas Equações ([[#eq-1|1]])-([[#eq-4|4]]), <math display="inline">f_m</math>, <math display="inline">g_m</math> e <math display="inline">S_m</math> também dependem destas últimas.  Por  outro lado, é preciso conhecer <math display="inline">[(q_m)_{xx}]_{i,j}^k</math> e <math display="inline">[(q_m)_{yy}]_{i,j}^k</math> previamente à obtenção de <math display="inline">(q_m)_{i,j}^k</math> no tempo <math display="inline">t_k</math>, conforme Eq. ([[#eq-10|10]]). A seção a seguir descreve a obtenção dessas  derivadas espaciais duplas.

====2.2.2 Avaliação de [(q<sub>m</sub>)<sub>xx</sub>]<sub>i,j</sub><sup>k</sup>, [(q<sub>m</sub>)<sub>xy</sub>]<sub>i,j</sub><sup>k</sup> e [(q<sub>m</sub>)<sub>yy</sub>]<sub>i,j</sub><sup>k</sup>====

Para prosseguir com o cálculo de <math display="inline">[(q_m)_{xx}]_{i,j}^k</math>, constrói-se primeiramente uma equação auxiliar a  partir da Eq. ([[#eq-5|5]]), derivando-a duas vezes em relação a <math display="inline">x</math>, obtendo as seguintes formas diferencial  e integral

<span id="eq-12"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\displaystyle \frac{\partial (q_m)_{xx}}{\partial{t}}+\frac{\partial }{\partial x}  \left[(f_m)_{xx}-(S_m)_x\right]+\frac{\partial (g_m)_{xx}}{\partial  y} = \textrm{div }\hat{\textbf{H}}_m = 0 \quad \textrm{e} \quad  \oint _{S(V)}\hat{\textbf{H}}_m \cdot d\textbf{s} = 0, \quad m = 1,2,3. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (12)
|}

As Equações em ([[#eq-12|12]]) possuem, conforme Equações ([[#eq-7|7]]) e ([[#eq-8|8]]), para todo <math display="inline">(x,y,t) \in </math> SE<math display="inline">(i,j,k)</math>, os seguintes análogos numéricos:

<span id="eq-13"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\textrm{div }\hat{\textbf{H}}_m^* = 0  \quad \textrm{e} \quad \oint _{S(\textrm{CE}(i,j,k))}\hat{\textbf{H}}_m^* \cdot  d\textbf{s} = 0, \quad \quad m = 1,2,3.   </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (13)
|}

O campo vetorial aproximado <math display="inline">\hat{\textbf{H}}_m^* = ((f_m^*)_{xx}-(S_m^*)_x,  (g_m^*)_{xx},(q_m^*)_{xx})</math> é constante, de modo que a avaliação da Eq.  ([[#eq-13|13]]) sobre o elemento de conservação descrito pela Fig. [[#img-1a|1a]] retorna, para <math display="inline">m = 1,2,3</math>:

<span id="eq-14"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>[(q_m)_{xx}]_{i,j}^k = \frac{1}{4}\left\{[(q_m)_{xx}]_{i-1/2,j}^{k-1/2}+[(q_m)_{xx}]_{i+1/2,j}^{k-1/2}+ [(q_m)_{xx}]_{i,j-1/2}^{k-1/2}+\right.</math><math> \left.[(q_m)_{xx}]_{i,j+1/2}^{k-1/2}\right\}+\frac{\Delta t}{2 \Delta x}\left\{ [(f_m)_{xx}]_{i-1/2,j}^{k-1/2}\right.</math>
|-
| style="text-align: center;" | <math>  \left.-[(f_m)_{xx}]_{i+1/2,j}^{k-1/2}  +[(S_m)_x]_{i+1/2,j}^{k-1/2}-[(S_m)_x]_{i-1/2,j}^{k-1/2}\right\}+</math><math>\frac{\Delta t}{2 \Delta y}\left\{ [(g_m)_{xx}]_{i,j-1/2}^{k-1/2}-[(g_m)_{xx}]_{i,j+1/2}^{k-1/2}\right\}.  </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (14)
|}

Para computar as variáveis de marcha <math display="inline">[(q_m)_{yy}]_{i,j}^k</math> e <math display="inline">[(q_m)_{xy}]_{i,j}^k</math> procede-se de modo  análogo, obtendo-se após todos os cálculos e simplificações:

<span id="eq-15"></span>
<span id="eq-16"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>{[}(q_m)_{yy}{]}_{i,j}^k = \frac{1}{4}\left\{[(q_m)_{yy}]_{i-1/2,j}^{k-1/2}+[(q_m)_{yy}]_{i+1/2,j}^{k-1/2}+ [(q_m)_{yy}]_{i,j-1/2}^{k-1/2}+\right.</math><math>\left.[(q_m)_{yy}]_{i,j+1/2}^{k-1/2}\right\}+\frac{\Delta t}{2 \Delta y}\left\{ [(g_m)_{yy}]_{i,j-1/2}^{k-1/2}\right.</math>
|-
| style="text-align: center;" | <math>   \left.-[(g_m)_{yy}]_{i,j+1/2}^{k-1/2}+[(S_m)_y]_{i,j+1/2}^{k-1/2}-[(S_m)_y] _ {  i,j-1/2} ^{k-1/2}\right\} +</math><math>\frac{\Delta t}{2 \Delta x}\left\{ [(f_m)_{yy}]_{i-1/2,j}^{k-1/2}-[(f_m)_{yy}]_{i+1/2,j}^{k-1/2}\right\}, </math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (15)
|-
| style="text-align: center;padding-top:15px" | <math> {[}(q_m)_{xy}{]}_{i,j}^k= \frac{1}{4}\left\{ [(q_m)_{xy}]_{i-1/2,j}^{k-1/2}+[(q_m)_{xy}]_{i+1/2,j}^{k-1/2}+ [(q_m)_{xy}]_{i,j-1/2}^{k-1/2}+[(q_m)_{xy}]_{i,j+1/2}^{k-1/2}\right\}</math>
|-
| style="text-align: center;" | <math>  +\frac{\Delta t}{2 \Delta x}\left\{ [(f_m)_{xy}]_{i-1/2,j}^{k-1/2}-[(f_m)_{xy}]_{i+1/2,j}^{k-1/2}+\frac{1}{2}[(S_m)_y]_{i+1/2,j}^{k-1/2}-\frac{1}{ 2}[(S_m)_y] _ { i-1/2 , j } ^{k-1/2}\right\}</math>
|-
| style="text-align: center;" | <math>   +\frac{\Delta t}{2 \Delta y}\left\{[(g_m)_{xy}]_{i,j-1/2}^{k-1/2}-[(g_m)_{xy}]_{i,j+1/2}^{k-1/2}  +\frac{1}{2}[(S_m)_x]_{i,j+1/2}^{k-1/2}-\frac{1}{2}[(S_m)_x] _ {i,j-1/2} ^{k-1/2}\right\}, \quad  m = 1,2,3.  </math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (16)
|}
|}

Aborda-se na seção a seguir o cálculo das derivadas de primeira ordem das variáveis dinâmicas <math display="inline">q_m</math>, <math display="inline">m =  1,2,3</math>.

====2.2.3 Avaliação de [(q<sub>m</sub>)<sub>x</sub>]<sub>i,j</sub><sup>k</sup> e [(q<sub>m</sub>)<sub>y</sub>]<sub>i,j</sub><sup>k</sup>====

As variáveis <math display="inline">[(q_m)_x]_{i,j}^k</math> e <math display="inline">[(q_m)_y]_{i,j}^k</math> podem ser determinadas por meio de uma estratégia  análoga àquela apresentada na seção [[#2.2.2 Avaliação de <math>[(q_m)_{xx}]_{i,j}^k</math>, <math>[(q_m)_{xy}]_{i,j}^k</math> e <math>[(q_m)_{yy}]_{i,j}^k</math>|2.2.2]] anterior, bastando, para isso, construir leis diferenciais  de conservação (já nas formas aproximadas)

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\displaystyle \frac{\partial (q_m^*)_{x}}{\partial{t}}+\frac{\partial (f_m^*)_{x}}{\partial x}  +\frac{\partial (g_m^*)_{x}}{\partial y} =\textrm{div }\check{\textbf{H}}_m^*= (S_m^*)_x, </math> <math>\quad  \frac{\partial (q_m^*)_{y}}{\partial{t}} +\frac{\partial (f_m^*)_{y}}{\partial x} +\frac{\partial  (g_m^*)_{y}}{\partial y} = \textrm{div }\tilde{\textbf{H}}_m^* =  (S_m^*)_y </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (17)
|}

e integrais

<span id="eq-18"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\oint _{S(\textrm{CE}(i,j,k))}\check{\textbf{H}}_m^*\cdot d\textbf{s}=  \int _{\textrm{CE}(i,j,k)}(S_m^*)_xdV \quad  \textrm{e} \quad  </math> <math>\oint _{S(\textrm{CE}(i,j,k))}\tilde{\textbf{H}}_m^*\cdot  d\textbf{s}=  \int _{\textrm{CE}(i,j,k)}(S_m^*)_ydV, \quad m = 1,2,3. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (18)
|}

As avaliações das equações em ([[#eq-18|18]]) fornecem as inclinações equacionadas por

<span id="eq-19"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>[(q_m)_{\xi }]_{i,j}^k = [Q_m^{(1)}(\Delta x,\xi )]_{i-1/2,j}^{k-1/2}+[Q_m^{(1)}(-\Delta  x,\xi )]_{i+1/2,j}^{k-1/2}+</math> <math>[Q_m^{(2)}(\Delta y,\xi )]_{i,j-1/2}^{k-1/2}+[Q_m^{(2)}(-\Delta  y,\xi )]_{i,j+1/2}^{k-1/2}, </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (19)
|}

para <math display="inline">\xi = x,y</math> e <math display="inline">m = 1,2,3</math>, em que

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>Q_m^{(1)}(\Delta x,\xi ) = \displaystyle  \frac{1}{4}\left\{(q_m)_{\xi }+\frac{\Delta  x}{4}(q_m)_{x\xi }+\frac{\Delta t}{4}\left[(f_m)_{x\xi }-(g_m)_{y\xi }\right]+\right.</math> <math>\left.\frac{2\Delta t}{\Delta  x}\left[(f_m)_{\xi }+\displaystyle \frac{\Delta t}{4}(f_m)_{t\xi }  \right]+\frac{\Delta t}{2} (S_m)_{\xi }\right\},</math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (20)
|-
| style="text-align: center;" | <math> Q_m^{(2)}(\Delta y,\xi ) = \displaystyle  \frac{1}{4}\left\{(q_m)_{\xi }+\frac{\Delta  y}{4}(q_m)_{y\xi }+\frac{\Delta t}{4}\left[(g_m)_{y\xi }-(f_m)_{x\xi }\right]+\right. </math> <math>\left.\frac{2\Delta t}{\Delta  y}\left[(g_m)_{\xi }+\displaystyle \frac{\Delta t}{4}(g_m)_{t\xi }  \right]+\frac{\Delta t}{2}(S_m)_{\xi }\right\}. </math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (21)
|}
|}

O esquema formado pelo conjunto de Equações ([[#eq-10|10]]),  ([[#eq-14|14]])-([[#eq-16|16]]) e  ([[#eq-19|19]]) é interessante apenas em problemas com soluções  de comportamento suave. Para problemas com descontinuidades ou formação de choque, é necessário modificar o  cálculo das variáveis em ([[#eq-14|14]])-([[#eq-16|16]]) e ([[#eq-19|19]]). A proposta apresentada  neste artigo é uma variação do procedimento adotado por <span id='citeF-38'></span><span id='citeF-39'></span>[[#cite-38|[38,39]]] e consiste num processo  iterativo.

Para prosseguir, considere as seguintes formas regressivas e progressivas:

<span id="eq-22"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>[(q_m^{\pm })_{\xi }]_{i,j}^k = \delta _{\xi }^{\pm }(q_m) \quad \textrm{e} \quad  [(q_m^{\pm })_{\zeta \xi }]_{i,j}^k = \delta _{\xi }^{\pm }[(q_m)_{\zeta }], \quad </math> <math> \xi = x, y, \quad \zeta = x,  y \quad \textrm{e} \quad m = 1,2,3 </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (22)
|}

tal que o operador <math display="inline">\delta _{\xi }^{\pm }</math> seja definido como

<span id="eq-23"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\delta _{\xi }^{\pm }(\phi ) =  \displaystyle \left\{\displaystyle \begin{array}{cl}\displaystyle \pm \frac{2}{\Delta \xi }\left[(\phi ')_{i\pm{1}/2,j}^k-\phi _{i,j}^k\right],&  \textrm{se } \;\xi = x\\ \displaystyle \pm \frac{2}{\Delta \xi }\left[(\phi ')_{i,j\pm{1}/2}^k-\phi _{i,j}^k\right],& \textrm{se }\;  \xi = y \end{array}\right.\quad \textrm{e} \quad </math> <math> \phi' = \left\{\begin{array}{cl}\displaystyle \phi +\frac{\Delta t}{2}\phi _{t}+\frac{\Delta t^2}{8}\phi _{tt}, & \textrm{se }\; \phi =  (q_m)_x, (q_m)_y, \\ \displaystyle \phi +\frac{\Delta t}{2}\phi _{t} , & \textrm{se }\; \phi =  (q_m)_{xx}, (q_m)_{xy}, (q_m)_{yy}. \end{array}\right. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (23)
|}

As aproximações para as derivadas serão definidas como

<span id="eq-24"></span>
<span id="eq-25"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>[(q_m)_{\xi }]_{i,j}^k = W\left\{[(q_m^+)_{\xi }]_{i,j}^k, [(q_m^-)_{\xi }]_{i,j}^k,  \alpha \right\}, </math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (24)
|-
| style="text-align: center;" | <math> {[}(q_m)_{\zeta \xi }{]}_{i,j}^k = W\left\{[(q_m^+)_{\zeta \xi }]_{i,j}^k,  [(q_m^-)_{\zeta \xi }]_{i,j}^k, \alpha \right\}, \quad </math> <math> m = 1,2,3, \quad \xi = x, y, \quad \zeta = x,  y \quad \textrm{e}  \quad  \alpha \geq 0,  </math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (25)
|}
|}

e a função <math display="inline">W</math> <span id='citeF-26'></span><span id='citeF-24'></span>[[#cite-26|[26,24]]]:

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>W(x_{-}, x_{+}, \alpha ) =   \left\{  \begin{array}{cl}\displaystyle \frac{|x_{+}|^{\alpha }x_{-}+|x_{-}|^{\alpha }x_{+}}{|x_{+}|^{\alpha }+|x_{-}|^{\alpha }}, &  \textrm{se} \quad |x_{+}|^{\alpha }+|x_{-}|^{\alpha } \neq 0, \\ 0, & \textrm{caso contrario.} \end{array} \right. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (26)
|}

Observe que neste caso, as Equações ([[#eq-24|24]]) e ([[#eq-25|25]]) representam médias ponderadas das diferenças  progressivas e regressivas descritas em ([[#eq-22|22]]) e ([[#eq-23|23]]). No caso de <math display="inline">\alpha = 0</math>,  ([[#eq-24|24]]) e ([[#eq-25|25]]) tornam-se diferenças finitas centrais <span id='citeF-23'></span>[[#cite-23|[23]]]. Adota-se neste trabalho <math display="inline">\alpha =1</math>.

Importante notar, no entanto, que as Equações ([[#eq-24|24]]) e ([[#eq-25|25]]) dependem da variável <math display="inline">(q_m)_{i,j}^k</math> que, por sua vez, depende de <math display="inline">[(q_m)_{xx}]_{i,j}^k</math> e <math display="inline">[(q_m)_{yy}]_{i,j}^k</math>, conforme Eq.  ([[#eq-10|10]]). Dito de outra forma, para obter <math display="inline">(q_m)_{i,j}^k</math> é necessário obter antes <math display="inline">[(q_m)_{xx}]_{i,j}^k</math> e <math display="inline">[(q_m)_{yy}]_{i,j}^k</math>, e vice-versa. Para sanar esta dificuldade, construímos o  processo iterativo descrito no Algoritmo [[#algorithm-1|1]].

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;font-size: 75%;">
'''Algoritmo 1:''' Algoritmo CE/SE para o cálculo das variáveis de marcha sobre um determinado ponto <math>(i,j,k)</math> da  malha.</div>
{| style="margin: 1em auto;border: 1px solid darkgray;font-size: 85%;"
|-
|
:'''início'''
::Defina uma tolerância <math>\tau </math>
::   Faça <math>m = 1</math>
::'''repita'''
:::Faça <math>l = 1</math>;
:::  Calcule <math>[(q_m)_{xx}]_{i,j}^k</math> e <math>[(q_m)_{yy}]_{i,j}^k</math> com as Equações ([[#eq-14|14]]) e  ([[#eq-15|15]]), respectivamente;
:::Calcule <math>(q_m^l)_{i,j}^k</math> com a  Eq. ([[#eq-10|10]]);
:::'''repita'''
::::Calcule <math>[(q_m)_x]_{i,j}^k</math> e <math>[(q_m)_y]_{i,j}^k</math> com a Eq. ([[#eq-24|24]]);
::::Calcule <math>[(q_m)_{xx}]_{i,j}^k</math>, <math>[(q_m)_{xy}]_{i,j}^k</math> e <math>[(q_m)_{yy}]_{i,j}^k</math> com  a Eq.([[#eq-25|25]]);
::::Calcule <math>(q_m^{l+1})_{i,j}^k</math> com a  Eq. ([[#eq-10|10]]);  
:::'''até''' <math>||(q_m^{l+1})_{i,j}^k-(q_m^l)_{i,j}^k||< \tau ;</math>      
::'''até''' <math>m = 3;</math>          
:'''fim'''
|
|}


A estabilidade do esquema CE/SE satisfaz a condição de Courant-Friedrichs-Lewy (CFL) que, para  as equações de Saint Venant unidimensional e bidimensional são as restrições <span id='citeF-23'></span>[[#cite-23|[23]]]

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>CFL = \displaystyle \frac{\Delta t}{\Delta x}\textrm{max}\left(|u|+\sqrt{gh}\right)\quad \textrm{e} \quad  </math> <math> CFL = \displaystyle \textrm{max}\left\{\frac{\Delta t}{\Delta x}\left(|u|+\sqrt{gh}\right),  \frac{\Delta t}{\Delta y} \left(|v|+\sqrt{gh}\right)\right\}, </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (27)
|}

respectivamente, onde o número de Courant satisfaz <math display="inline">0 < CFL \leq 1</math>. Maiores detalhes sobre  estabilidade do método CE/SE pode ser encontrado em <span id='citeF-39'></span>[[#cite-40|[40]]].

==3 Exemplos Numéricos==

Esta seção tem por objetivo avaliar e validar o esquema numérico desenvolvido. Para isso, alguns problemas  teste são avaliados, como testes de acurácia em uma e duas dimensões, bem como problemas clássicos de ruptura  de barragem unidimensionais (subseção [[#3.1 Exemplos Unidimensionais|3.1]]) e bidimensionais (subseção [[#3.2 Exemplos Bidimensionais|3.2]]).

===3.1 Exemplos Unidimensionais===

<span id='theorem-ex1'></span>'''Exemplo 3.1.1''':  O objetivo deste exemplo teste é verificar experimentalmente a ordem de acurácia do esquema  CE/SE. Considere o sistema hiperbólico linear unidimensional

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>  \textbf{U}_t+\textbf{F}(\textbf{U})_x = \textbf{0}, \quad  \textrm{onde} \quad \textbf{U} = (\phi , u), \quad \textbf{F} = (cu, c\phi ) </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (28)
|}

e <math>c</math> é uma constante dada. Impondo as condições inicias <math>\phi (x,0) =  -c^{-1}\sin (2\pi x)</math> e <math>u(x,0) = 0</math>, é possível obter a seguinte solução exata <span id='citeF-40'></span><span id='citeF-41'></span>[[#cite-41|[41,42]]]:

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>  \phi (x,t) = -c^{-1}\cos (2\pi c t)\sin (2\pi x) \quad \textrm{e} \quad </math>  <math> u(x,t) = c^{-1}\sin (2\pi ct)\cos (2\pi  x). </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (29)
|}

As soluções são computadas com as Equações ([[#eq-10|10]]), ([[#eq-14|14]])-([[#eq-16|16]]) e  ([[#eq-19|19]]). A  Tabela [[#table-1|1]] apresenta os erros calculados nas normas <math>L^1</math> e <math>L^2</math>, nas simulações  realizadas sobre o domínio computacional <math>-1 \leq x \leq 1</math>, <math>c = 1</math>, <math>CFL = 0.8</math>  até atingir o tempo <math>t = 1.2</math>s, sendo <math>n</math> o número de células no espaço. O experimento verifica a terceira  ordem de acurácia do esquema para ambas as variáveis <math>u</math> e <math>\phi </math>. 


<div class="center" style="font-size: 75%;">'''Tabela 1'''. Ordem de acurácia experimental do método CE/SE para o sistema hiperbólico linear  unidimensional.</div>

{| style="width: 100%;margin: 1em auto 0.1em auto;border-collapse: collapse;font-size:85%" 
|- style="border-top: 1px solid;"  
|
| colspan='4' style=" text-align: center;border-bottom: 1px solid;" | Norma <math>L^1</math> 
| colspan='1' style=" width: 20px;" |
| colspan='4' style=" text-align: center;border-bottom: 1px solid;" |  Norma <math>L^2</math>
|- 
| style="text-align: center;"| ''n''
|  style="text-align: center;"|Erro (<math>\phi </math>) 
| style="text-align: center;"|Ordem 
|  style="text-align: center;"|Erro (<math>u</math>)
|  style="text-align: center;"|Ordem  
| 
|  style="text-align: center;"|Erro (<math>\phi </math>)
| style="text-align: center;"|Ordem 
|  style="text-align: center;"|Erro  (<math>u</math>)
|  style="text-align: center;"| Ordem 
|- style="border-top: 1px solid;"
| style="text-align: center;"| 20 
| style="text-align: center;"|8.322<math>\times 10^{-3}</math> 
| style="text-align: center;"|&#8211; 
| style="text-align: center;"|1.662<math>\times 10^{-2}</math> 
| style="text-align: center;"|&#8211; 
|  
| style="text-align: center;"|7.157<math>\times 10^{-3}</math> 
|style="text-align: center;"| &#8211; 
| style="text-align: center;"|1.347<math>\times 10^{-2}</math> 
| style="text-align: center;"|&#8211; 
|-
|  style="text-align: center;"|40 
| style="text-align: center;"|1.849<math>\times 10^{-3}</math> 
| style="text-align: center;"| 2.170 
| style="text-align: center;"|2.524<math>\times 10^{-3}</math> 
| style="text-align: center;"|2.719
|  
| style="text-align: center;"|1.488<math>\times 10^{-3}</math>  
| style="text-align: center;"|2.266 
| style="text-align: center;"|2.042<math>\times 10^{-3}</math> 
| style="text-align: center;"|2.722
|-
| style="text-align: center;"| 80 
| style="text-align: center;"|3.067<math>\times 10^{-4}</math> 
| style="text-align: center;"| 2.592 
| style="text-align: center;"|3.401<math>\times 10^{-4}</math> 
|style="text-align: center;"| 2.892
|  
| style="text-align: center;"|2.443<math>\times 10^{-4}</math>  
| style="text-align: center;"|2.607 
| style="text-align: center;"|2.715<math>\times 10^{-4}</math> 
| style="text-align: center;"|2.911
|-
|  style="text-align: center;"|160 
| style="text-align: center;"|4.369<math>\times 10^{-5}</math> 
|  style="text-align: center;"|2.811 
| style="text-align: center;"|4.386<math>\times 10^{-5}</math> 
| style="text-align: center;"|2.955
|  
| style="text-align: center;"|3.493<math>\times 10^{-5}</math>  
| style="text-align: center;"|2.806 
| style="text-align: center;"|3.465<math>\times 10^{-5}</math> 
| style="text-align: center;"|2.970
|-
|  style="text-align: center;"|320 
| style="text-align: center;"|5.813<math>\times 10^{-6}</math> 
|  style="text-align: center;"|2.910 
| style="text-align: center;"|5.566<math>\times 10^{-6}</math> 
| style="text-align: center;"|2.978
|  
| style="text-align: center;"|4.665<math>\times 10^{-6}</math>  
| style="text-align: center;"|2.904 
| style="text-align: center;"|4.369<math>\times 10^{-6}</math> 
| style="text-align: center;"|2.987
|- style="border-bottom: 1px solid;"
| style="text-align: center;"| 640 
| style="text-align: center;"|7.490<math>\times 10^{-7}</math> 
| style="text-align: center;"| 2.956 
| style="text-align: center;"|7.010<math>\times 10^{-7}</math> 
| style="text-align: center;"|2.989
|  
| style="text-align: center;"|6.030<math>\times 10^{-7}</math>  
| style="text-align: center;"|2.952 
|style="text-align: center;"| 5.490<math>\times 10^{-7}</math> 
| style="text-align: center;"|2.992
|}


<span id='theorem-fig_DBS01'></span>'''Exemplo 3.1.2''':  Considera-se agora um problema ideal de ruptura de barragem sobre um domínio molhado, isto é, a quebra  de barragem é instantânea, o fundo é plano e não existe resistência ao escoamento. As condições iniciais para  esta configuração seguem o clássico problema de Riemann

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>  h(x,0) =   \left\{   \begin{array}{cl}    h_l, & \textrm{ para } 0 \textrm{m} \leq x \leq x_0 \textrm{m}\\    h_r , & \textrm{ para } x_0 \textrm{m} < x \leq L \textrm{m}   \end{array}  \right., \quad </math> <math>\textrm{ com } \quad h_l \geq h_r \quad \textrm{e} \quad u(x,0) = 0 \textrm{m/s}. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (30)
|}

O domínio considerado é <math>0\mbox{m}</math> <math>\leq x \leq 2000\mbox{m}</math>, <math>x_0 = 1000\mbox{m}</math> e a solução analítica  pode ser encontrada em <span id='citeF-42'></span>[[#cite-43|[43]]] ou <span id='citeF-43'></span>[[#cite-44|[44]]]. Utiliza-se, para fins de análise  de comportamento, a mesma relação de equações do Ex. [[#theorem-ex1|3.1.1]] anterior, isto é, Equações ([[#eq-3.10|10]]),  ([[#eq-14|14]])-([[#eq-16|16]]) e ([[#eq-19|19]]). As Fig. [[#theorem-fig_DBS01|2a]] e  Fig. [[#img-2a|2b]] demonstram o comportamento da solução numérica em relação à analítica, calculadas no  tempo <math>t = 52\mbox{s}</math>, com <math>n+1 = 201</math> pontos, incremento espacial <math>\Delta x = 10\mbox{m}</math>, <math>CFL = 0.8</math> e  profundidades iniciais a montante e a jusante <math>h_l = 10\mbox{m}</math> e <math>h_r = 5\mbox{m}</math>, respectivamente. Note  que a resposta numérica é coerente com a analítica, embora apresente suavidade que a distancie nas regiões  com mudanças abruptas.   


<div id='img-2a'></div>
<div id='img-2b'></div>
{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|style="text-align: center;padding-top:10px"|[[File:Draft_Aparicio_Nogué_501194349_7076_fig2a.png]]
|style="text-align: center;padding-top:10px"|
[[File:Draft_Aparicio_Nogué_501194349_1886_fig2b.png]]
|-
|style="text-align: center;font-size: 75%;"| (a)
|style="text-align: center;font-size: 75%;"| (b)
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding-top:10px;" | '''Figura 2'''. Elevação da superfície da água <math>h</math> e velocidade <math>u</math> para o problema de ruptura de  barragem, computada no tempo <math>t = 52\mbox{s}</math>, utilizando uma malha uniforme com <math>201</math>  pontos, incremento espacial <math>\Delta x = 10\mbox{m}</math> e <math>h_r/h_l = 0.5</math>.
|}


Utilizando o Algoritmo [[#algorithm-1|1]], sob as mesmas condições e com os mesmos parâmetros, obtém-se os  gráficos constantes na  Fig. [[#img-3|3]]. Observa-se que esta solução numérica é superior à anterior,  sobretudo no que tange às regiões de rápidas mudanças.  Conforme Zhang ''et al.'' (2012) <span id='citeF-23'></span>[[#cite-23|[23]]],  a razão <math>h_r/h_l</math> é um importante índice para  julgar a aplicabilidade e a acurácia de esquemas numéricos no modelo 1D de ruptura de barragem. Segundo os  mesmos autores, os regimes de escoamentos subcrítico e supercrítico existem simultaneamente num canal  sem fricção, horizontal e retangular, quando <math>h_r/h_l < 0.138</math>.  Altera-se, neste sentido, estes  parâmetros para uma razão <math>h_r/h_l = 0.001</math>. Os resultados simulados são mostrados nas Figs. [[#img-3c|3c]] e [[#img-3a|3d]].   


<div id='img-3a'></div>
<div id='img-3b'></div>
<div id='img-3c'></div>
<div id='img-3d'></div>
{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|style="text-align: center;padding-top:10px"|
[[File:Draft_Aparicio_Nogué_501194349_4710_fig3a.png]]
|style="text-align: center;padding-top:10px"|
[[File:Draft_Aparicio_Nogué_501194349_9130_fig3b.png]]
|-
|style="text-align: center;font-size: 75%;"| (a)
|style="text-align: center;font-size: 75%;"| (b)
|- 
|style="text-align: center;padding-top:10px"|
[[File:Draft_Aparicio_Nogué_501194349_7053_fig3c.png]]
|style="text-align: center;padding-top:10px"|
[[File:Draft_Aparicio_Nogué_501194349_1821_fig3d.png]]
|-
|style="text-align: center;font-size: 75%;"| (c)
|style="text-align: center;font-size: 75%;"| (d)
|- 
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding-top:10px;"| '''Figura 3'''. Elevação da superfície da água <math>h</math> e velocidade <math>u</math> para o problema de ruptura de  barragem, no tempo <math>t = 52\mbox{s}</math>, numa malha uniforme com <math>201</math> pontos e respectivo incremento espacial <math>\Delta x = 10\mbox{m}</math>, calculados com o Algoritmo [[#algorithm-1|1]], sendo: (a)-(b) <math>h_l = 10\mbox{m}</math> e <math>h_r = 5\mbox{m}</math>, (c)-(d) <math>h_l = 10\mbox{m}</math> e <math>h_r = 0.01\mbox{m}</math>.
|}

===3.2 Exemplos Bidimensionais===

<span id='theorem-tabela02'></span>'''Exemplo 3.2.1''':  Considere agora o sistema hiperbólico linear bidimensional

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>  \textbf{U}_t+\textbf{F}(\textbf{U})_x  +\textbf{G}(\textbf{U})_y= \textbf{0}, \quad  \textrm{onde} \quad \textbf{U} = (\phi , u,v), </math><math> \quad \textbf{F} = (cu, c\phi ,0), \quad  \textbf{G} = (cv, 0,c\phi ) </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (31)
|}

e <math>c</math> é uma constante dada. Para os dados iniciais <math>\phi (x,y,0) =  -c^{-1}[\sin (2\pi x)+\sin (2\pi y)]</math> e <math>u(x,y,0) = v(x,y,0) = 0</math>, este problema admite a seguinte solução  exata <span id='citeF-40'></span><span id='citeF-41'></span>[[#cite-40|[41,42]]]:

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>  \phi = -c^{-1}\cos (2\pi c t)[\sin (2\pi x)+\sin (2\pi y)], \quad </math><math> u = c^{-1}\sin (2\pi  ct)\cos (2\pi x) \quad \textrm{e} </math><math> \quad  v = c^{-1}\sin (2\pi ct)\cos (2\pi y). </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (32)
|}


A Tabela [[#theorem-tabela02|2]] apresenta os erros numéricos calculados nas normas <math>L^1</math> e <math>L^2</math>, no tempo <math>t =  0.2</math>s, com parâmetros especificados em <math>c = 1</math>, <math>CFL = 0.4</math> e domínio computacional definido em <math>-1 \leq x, y  \leq 1</math>, sendo <math>n_x \times n_y</math> o número de células utilizadas na discretização. O experimento confirma  novamente a acurácia de terceira ordem do esquema CE/SE.


<div class="center" style="font-size: 75%;">'''Tabela 2'''. Ordem de acurácia experimental do método CE/SE para o sistema hiperbólico linear  bidimensional..</div>

{| style="width: 100%;margin: 1em auto 0.1em auto;border-collapse: collapse;font-size:85%" 
|- style="border-top: 1px solid;"  
|
| colspan='4' style=" text-align: center;border-bottom: 1px solid;" | Norma <math>L^1</math> 
| colspan='1' style=" width: 20px;" |
| colspan='4' style=" text-align: center;border-bottom: 1px solid;" |  Norma <math>L^2</math>
|- 
| style="text-align: center;"| <math>n_x \times n_y</math> 
| style="text-align: center;"| Erro (<math>\phi </math>) 
| style="text-align: center;"|Ordem 
| style="text-align: center;"| Erro (<math>u, v</math>)
| style="text-align: center;"| Ordem  
| 
| style="text-align: center;"| Erro (<math>\phi </math>)
| style="text-align: center;"|Ordem 
| style="text-align: center;"| Erro  (<math>u, v</math>)
|  style="text-align: center;"| Ordem 
|- style="border-top: 1px solid;"
| style="text-align: center;"|<math>20\times 20</math> 
| style="text-align: center;"| 2.595<math>\times 10^{-3}</math> 
| style="text-align: center;"|&#8211; 
| style="text-align: center;"| 1.283<math>\times 10^{-2}</math> 
| style="text-align: center;"|&#8211; 
|  
| style="text-align: center;"|1.797<math>\times 10^{-3}</math> 
| style="text-align: center;"|&#8211; 
|  style="text-align: center;"|7.478<math>\times 10^{-3}</math> 
| style="text-align: center;"|&#8211;  
|-
|style="text-align: center;"| <math>40\times 40</math> 
|  style="text-align: center;"|4.231<math>\times 10^{-4}</math> 
| style="text-align: center;"| 2.617 
| style="text-align: center;"| 2.032<math>\times 10^{-3}</math> 
|  style="text-align: center;"|2.658 
|  
| style="text-align: center;"|3.377<math>\times 10^{-4}</math> 
| style="text-align: center;"|  2.412 
| style="text-align: center;"| 1.169<math>\times 10^{-3}</math> 
| style="text-align: center;"| 2.678 
|-
| style="text-align: center;"|<math>80\times 80</math> 
| style="text-align: center;"| 7.519<math>\times 10^{-5}</math> 
| style="text-align: center;"| 2.492 
|  style="text-align: center;"|2.924<math>\times 10^{-4}</math> 
| style="text-align: center;"| 2.797 
|  
| style="text-align: center;"|6.111<math>\times 10^{-5}</math> 
|  style="text-align: center;"| 2.466 
| style="text-align: center;"| 1.653<math>\times 10^{-4}</math> 
|  style="text-align: center;"|  2.822 
|-
| style="text-align: center;"|<math>160 \times 160</math> 
| style="text-align: center;"| 1.197<math>\times 10^{-5}</math> 
| style="text-align: center;"| 2.651 
|  style="text-align: center;"|3.938<math>\times 10^{-5}</math> 
| style="text-align: center;"| 2.892 
|  
| style="text-align: center;"|9.602<math>\times  10^{-6}</math> 
| style="text-align: center;"| 2.670 
| style="text-align: center;"| 2.217<math>\times 10^{-5}</math> 
|  style="text-align: center;"|2.898 
|-
| style="text-align: center;"|<math>320\times 320</math> 
|  style="text-align: center;"|1.735<math>\times 10^{-6}</math> 
|  style="text-align: center;"|2.786 
|  style="text-align: center;"|5.128<math>\times 10^{-6}</math> 
|  style="text-align: center;"|2.941 
|  
| style="text-align: center;"|1.397<math>\times 10^{-6}</math> 
| style="text-align: center;"|  2.781 
|  style="text-align: center;"|2.895<math>\times 10^{-6}</math> 
| style="text-align: center;"| 2.937 
|- style="border-bottom: 1px solid;"
| style="text-align: center;"|<math>640\times 640</math> 
| style="text-align: center;"| 2.364<math>\times 10^{-7}</math> 
| style="text-align: center;"| 2.876 
| style="text-align: center;"| 6.554<math>\times 10^{-7}</math> 
|  style="text-align: center;"|2.968 
|  
| style="text-align: center;"|1.939<math>\times 10^{-7}</math> 
|  style="text-align: center;"| 2.849 
|  style="text-align: center;"|3.723<math>\times 10^{-7}</math> 
| style="text-align: center;"| 2.959 
|}



<span id='theorem-fig_RB2D'></span>'''Exemplo 3.2.2''': Este problema hipotético bidimensional é um exemplo  utilizado por <span id='citeF-23'></span><span id='citeF-44'></span><span id='citeF-45'></span>[[#cite-23|[23,45,46]]]. Neste problema as velocidades iniciais são todas  nulas, a profundidade a montante é de <math display="inline">10</math>m, enquanto que a profundidade a jusante é  assumida ser <math display="inline">5</math>m ou <math display="inline">0.1</math>m. O domínio computacional consiste de uma região de <math display="inline">200\mbox{m}\times  200\mbox{m}</math>, com uma parede que se estende paralelamente ao eixo <math display="inline">y</math>, tendo <math display="inline">10</math>m de largura e está centrada  em <math display="inline">x = 100</math>m. A falha é suposta ser instantânea, possui <math display="inline">75</math>m de extensão a partir de <math display="inline">y = 95</math>m. O canal é  horizontal e desconsidera-se resistência ao escoamento. Espera-se a formação de uma frente  de choque após o rompimento. Os resultados em <math display="inline">t = 7.2</math>s são computados a  partir de uma malha uniforme composta por <math display="inline">101 \times 101</math> células e <math display="inline">CFL = 0.4</math>. Constam nas Figuras  [[#img-4a|4a]] e [[#img-4b|4b]] gráficos da profundidade, vetores velocidade e curvas de nível para o  problema com profundidade inicial a jusante de <math display="inline">5</math>m, enquanto que as Figuras [[#img-4c|4c]] e  [[#img-4d|4d]] referem-se ao problema com profundidade inicial a jusante de <math display="inline">0.1</math>m. Os resultados são  consistentes com aqueles presentes na literatura <span id='citeF-23'></span><span id='citeF-44'></span><span id='citeF-45'></span>[[#cite-23|[23,45,46]]].

<div id='img-4a'></div>
<div id='img-4b'></div>
<div id='img-4c'></div>
<div id='img-4d'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|style="text-align: center;" |[[Image:draft_Aparicio Nogué_501194349-fig_RB2D.png|570px|]]
|style="text-align: center;" |[[Image:draft_Aparicio Nogué_501194349-fig_CV_RB2D.png|400px|]]
|-
|style="text-align: center;font-size:75%;" |(a)
|style="text-align: center;font-size:75%;" |(b)
|-
|style="text-align: center;" |[[Image:draft_Aparicio Nogué_501194349-fig_RB2D1.png|570px|]]
|style="text-align: center;" |[[Image:draft_Aparicio Nogué_501194349-fig_CV1_RB2D.png|400px|]]
|-
|style="text-align: center;font-size:75%;" |(c)
|style="text-align: center;font-size:75%;" |(d)
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding-top:10px;"  | '''Figura 4'''. Elevação da superfície da água <math>h</math>, (a) e (c), vetores velocidade e curvas de nível, (b) e  (d), para a solução calculada no tempo <math>t = 7.2\mbox{s}</math> com <math>CFL = 0.4</math>, do problema de ruptura de  barragem anti-simétrica em um domínio horizontal e sem fricção, com profundidades a montante de <math>10</math>m e a  jusante de <math>5</math>m (a)-(b) e <math>0.1</math>m (c)-(d).
|}

<span id='theorem-fig_RBC2D'></span>'''Exemplo 3.2.3''': Este problema teste visa avaliar a  habilidade do esquema em preservar simetria. Considera-se um domínio de <math display="inline">50\mbox{m}\times 50\mbox{m}</math> com  condições iniciais:

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>h(x,y,0) = \left\{   \begin{array}{cl}10\mbox{m}, & \mbox{ se } (x-25)^2+(y-25)^2\leq 10.5^2\\     1\mbox{ m}, & \mbox{ caso contrario}   \end{array} \right.\quad </math><math>\textrm{e} \quad \quad u(x,y,0) = v(x,y,0) = 0. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (33)
|}

No instante da falha da barragem, supõe-se que a parede circular seja removida completamente e,  subsequentemente, formam-se ondas que se espalham radialmente. A solução numérica é computada numa malha  retangular uniforme com <math display="inline">101\times 101</math> células e o passo de tempo é tal que <math display="inline">CFL = 0.6</math>. A Fig.  [[#img-5a|5a]] apresenta o perfil da superfície da água <math display="inline">0.72</math>s após a hipotética falha na barragem  circular. Vetores velocidade e curvas nível relativas à superfície <math display="inline">h</math> estão dispostos na Fig.  [[#img-5b|5b]]. A simetria da solução numérica é bem preservada e está de acordo com a  literatura <span id='citeF-7'></span><span id='citeF-46'></span>[[#cite-7|[7,47]]].

<div id='img-5a'></div>
<div id='img-5b'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|style="text-align: center;" |[[Image:draft_Aparicio Nogué_501194349-fig_RBC2D.png|570px|]]
|style="text-align: center;" |[[Image:draft_Aparicio Nogué_501194349-fig_CV1_RBC2D.png|400px|]]
|-
|style="text-align: center;font-size:75%;" |(a)
|style="text-align: center;font-size:75%;" |(b)
|-
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding-top:10px;" | '''Figura 5'''. Perfil da superfície da água <math>h</math> (a), vetores velocidade e curvas de nível  (b), para a  solução calculada no tempo <math>t = 0.72\mbox{s}</math> com <math>CFL = 0.6</math>, do problema de ruptura de  barragem circular em um domínio horizontal e sem fricção.
|}

==4. Conclusões==

Este artigo apresentou o desenvolvimento de um novo esquema CE/SE explícito para a solução das equações de  águas rasas em uma e duas dimensões. Na formulação construída, as variáveis dinâmicas e, consequentemente, as  leis diferencial e integral de conservação, foram aproximadas localmente por expansões de Taylor de segunda  ordem. O conjunto de variáveis de marcha tornou-se constituído pelas variáveis fluxo <math display="inline">q_m</math> e suas  derivadas espaciais de primeira <math display="inline">(q_m)_x</math>, <math display="inline">(q_m)_y</math> e segunda ordem <math display="inline">(q_m)_{xx}, (q_m)_{xy}</math> e <math display="inline">(q_m)_{yy}</math>.  Os experimentos computacionais realizados indicaram acurácia de terceira ordem sobre os problemas  hiperbólicos lineares testados. As soluções numéricas dos problemas teste de ruptura de barragem  unidimensional e bidimensional, caracterizados pela formação de choque, descontinuidade ou simetria,  apresentaram concordância com as soluções analíticas e/ou com a literatura. Conclui-se, com isso, que o  esquema proposto possui considerável habilidade em capturar choques e descontinuidades fazendo com que seja  uma boa ferramenta para análise de fluxo de ruptura de barragem.

==Agradecimentos==

Os autores agradecem à CAPES pelo apoio financeiro, ao Programa de Pós-Graduação em Métodos Numéricos em Engenharia (PPGMNE - UFPR) e ao Instituto Federal de Educação, Ciência e Tecnologia Catarinense, Campus Araquari, pelo apoio à pesquisa.

==Referências==

<div class="auto" style="width: auto; margin-left: auto; margin-right: auto;font-size: 85%;">
<div id=cite-1></div>
[[#citeF-1|[1]]] S. Yang, A. Kurganov, Y. Liu, Well-balanced central schemes on overlapping cells with constant subtraction techniques for the Saint-Venant
shallow water system. Journal of Scientific Computing 63 (3) (2015) 678–698.

<div id=cite-2></div>
[[#citeF-2|[2]]] J. Singh, M. S. Altinakar, Y. Ding, Two-dimensional numerical modeling of dam-break flows over natural terrain using a central explicit
scheme. Advances in Water Resources 34 (10) (2011) 1366–1375.

<div id=cite-3></div>
[[#citeF-3|[3]]] J. Soler, E. Bladé, M. Sánchez-Juny, Ensayo comparativo entre modelos unidimensionales y bidimensionales en la modelización de la rotura
de una balsa de materiales sueltos erosionables. Revista Internacional de Métodos Numéricos para Cálculo y Diseño en Ingeniería 28 (2)
(2012) 103–111.

<div id=cite-4></div>
[[#citeF-4|[4]]] H. Ozmen-Cagatay, S. Kocaman, H. Guzel, Investigation of dam-break flood waves in a dry channel with a hump. Journal of Hydroenvironment
Research 8 (3) (2014) 304–315.

<div id=cite-5></div>
[[#citeF-5|[5]]] S. Kocaman, H. Ozmen-Cagatay, Investigation of dam-break induced shock waves impact on a vertical wall, Journal of Hydrology 525 (2015)
1–12.

<div id=cite-6></div>
[[#citeF-6|[6]]] T.-J. Chang, H.-M. Kao, K.-H. Chang, M.-H. Hsu, Numerical simulation of shallow-water dam break flows in open channels using smoothed
particle hydrodynamics. Journal of Hydrology 408 (1-2) (2011) 78–90.

<div id=cite-7></div>
[[#citeF-7|[7]]] R. Akoh, S. Li, F. Xiao, A CIP/multi-moment finite volume method for shallow water equations with source terms. Int. J. Numer. Meth.
Fluids 56 (2008) 2245–2270.

<div id=cite-8></div>
[[#citeF-8|[8]]] K. Toda, Y. Ogata, T. Yabe, Multi-dimensional conservative semi-Lagrangian method of characteristics CIP for the shallow water equations.
Journal of Computational Physics 228 (13) (2009) 4917–4944.

<div id=cite-9></div>
[[#citeF-9|[9]]] R. Borsche, A well-balanced solver for the Saint Venant equations with variable cross-section. Journal of Numerical Mathematics 23 (2)
(2015) 99–115.

<div id=cite-10></div>
[[#citeF-10|[10]]] J.-W.Wang, R.-X. Liu, Combined finite volume-finite element method for shallow water equations. Computers & Fluids 34 (10) (2005) 1199–1222.

<div id=cite-11></div>
[[#citeF-11|[11]]] M. C. Díaz, J. López-García, C. Parés, High order exactly well-balanced numerical methods for shallow water systems. Journal of Computational Physics 246 (2013) 242–264.

<div id=cite-12></div>
[[#citeF-12|[12]]] F. Benkhaldoun, S. Sari, M. Seaid, Projection finite volume method for shallow water flows. Mathematics and Computers in Simulation 118
(2015) 87 – 101, MAMERN V-2013: 5th International Conference on Approximation Methods and Numerical Modeling in Environment and Natural Resources.

<div id=cite-13></div>
[[#citeF-13|[13]]] S.-J. Liang, T.-W. Hsu, Least-squares finite-element method for shallow-water equations with source terms. Acta Mechanica Sinica 25 (5)
(2009) 597–610.

<div id=cite-14></div>
[[#citeF-14|[14]]] C. Dawson, J. H. Videman, A streamline diffusion finite element method for the viscous shallow water equations. Journal of Computational
and Applied Mathematics 251 (2013) 1–7.

<div id=cite-15></div>
[[#citeF-15|[15]]] G. Kesserwani, D. Caviedes-Voullième, N. Gerhard, S. M Multiwavelet discontinuous Galerkin h-adaptive shallow water model. Computer
Methods in Applied Mechanics and Engineering 294 (2015) 56–71.

<div id=cite-16></div>
[[#citeF-16|[16]]] D. Wirasaet, S. Brus, C. Michoski, E. Kubatko, J. Westerink, C. Dawson, Artificial boundary layers in discontinuous Galerkin solutions to
shallow water equations in channels. Journal of Computational Physics 299 (2015) 597–612.

<div id=cite-17></div>
[[#citeF-17|[17]]] M. Rasulov, Z. Aslan, O. Pakdil, Finite differences method for shallow water equations in a class of discontinuous functions. Applied
Mathematics and Computation 160 (2) (2005) 343–353.

<div id=cite-18></div>
[[#citeF-18|[18]]] X. Zhengfu, S. Chi-Wang, Anti-diffusive finite difference weno methods for shallow water with transport of pollutant. Journal of Computational Mathematics 24 (3) (2006) 239–251.

<div id=cite-19></div>
[[#citeF-19|[19]]] Y.-H. Hwang, A characteristic particle method for the Saint Venant equations. Computers & Fluids 76 (2013) 58–72.

<div id=cite-20></div>
[[#citeF-20|[20]]] H. Liu, H. Wang, S. Liu, C. Hu, Y. Ding, J. Zhang, Lattice Boltzmann method for the Saint-Venant equations. Journal of Hydrology 524
(2015) 411–416.

<div id=cite-21></div>
[[#citeF-21|[21]]] S.-C. Chang, W.-W. To, A new numerical framework for solving conservation laws - the method of space-time conservation element and
solution element. Technical Memo TM 104495, NASA, Lewis Research Center, NASA (August 1991).

<div id=cite-22></div>
[[#citeF-22|[22]]] T. Molls, F. Molls, Space-time conservation method applied to Saint Venant equations. Journal of Hydraulic Engineering 124 (5) (1998)
501–508.

<div id=cite-23></div>
[[#citeF-23|[23]]] Y. Zhang, Z. Zeng, J. Chen, The improved space-time conservation element and solution element scheme for two-dimensional dam-break
flow simulation. International Journal for Numerical Methods in Fluids 68 (5) (2012) 605–624.

<div id=cite-24></div>
[[#citeF-24|[24]]] S. Qamar, S. Zia, W. Ashraf, The space-time CE/SE method for solving single and two-phase shallow flow models. Computers & Fluids 96
(2014) 136–151.

<div id=cite-25></div>
[[#citeF-25|[25]]] S.-C. Chang, New developments in the method of space-time conservation element and solution element - aplications to the Euler and
Navier-Stokes equations. Technical Memo TM 106226, NASA, Lewis Research Center, NASA (August 1993).

<div id=cite-26></div>
[[#citeF-26|[26]]] S.-C. Chang, The method of space-time conservation element and solution elemen - a new approach for solving the Navier-Stokes and Euler
equations. Journal of Computational Physics 119 (2) (1995) 295–324.

<div id=cite-27></div>
[[#citeF-27|[27]]] S.-C. Chang, X.-Y. Wang, C.-Y. Chow, The space-time conservation element and solution element method: A new high-resolution and
genuinely multidimensional paradigm for solving conservation laws. Journal of Computational Physics 156 (1) (1999) 89–136.

<div id=cite-28></div>
[[#citeF-28|[28]]] S.-C. Chang, X.-Y. Wang, W.-M. To, Application of the space-time conservation element and solution element method to one-dimensional
convection-diffusion problems. Journal of Computational Physics 165 (1) (2000) 189–215.

<div id=cite-29></div>
[[#citeF-29|[29]]] Z.-C. Zhang, S. T. John Yu, X.-Y.Wang, S.-C. Chang, A. Himansu, P. C. E. Jorgenson, The CE/SE method for Navier-Stokes equations using
unstructured meshes for flows at all speeds. AIAA 2000-0393 (2000) 7.

<div id=cite-30></div>
[[#citeF-30|[30]]] Z.-C. Zhang, S. T. J. Yu, S.-C. Chang, A space-time conservation element and solution element method for solving the two- and threedimensional unsteady Euler equations using quadrilateral and hexahedral meshes. J. Comput. Phys. 175 (1) (2002) 168–199.

<div id=cite-31></div>
[[#citeF-31|[31]]] S.-C. Chang, A new approach for constructing highly stable high order CESE schemes. Technical Memo TM 2010-216766, NASA, Lewis
Research Center, NASA (August 2010).

<div id=cite-32></div>
[[#citeF-32|[32]]] D. L. Bilyeu, S.-T. J. Yu, Y. Y. Chen, J. L. Cambier, A two-dimensional fourth-order unstrutured-meshed Euler solver based on the CESE
method. J. Comput. Phys. 257 (1) (2014) 981–999.

<div id=cite-33></div>
[[#citeF-33|[33]]] S. Jerez, J. Romero, M. Roselló, F. Arnau, A semi-implicit space-time CE-SE method to improve mass conservation through tapered ducts in
internal combustion engines. Mathematical and Computer Modelling 40 (9-10) (2004) 941–951.

<div id=cite-34></div>
[[#citeF-34|[34]]] S. Qamar, G. Warnecke, Application of space-time CE/SE method to shallow water magnetohydrodynamic equations. Journal of Computational
and Applied Mathematics 196 (1) (2006) 132–149.

<div id=cite-35></div>
[[#citeF-35|[35]]] X. Luo, M. Wang, J. Yang, G. Wang, The space-time CESE method applied to phase transition of water vapor in compressible flows.
Computers & Fluids 36 (7) (2007) 1247–1258.

<div id=cite-36></div>
[[#citeF-36|[36]]] Y. Chou, R.-J. Yang, Application of CESE method to simulate non-Fourier heat conduction in finite medium with pulse surface heating.
International Journal of Heat and Mass Transfer 51 (13-14) (2008) 3525–3534.

<div id=cite-37></div>
[[#citeF-37|[37]]] Z.-C. Zhang, S. T. John Yu, S.-C. Chang, A. Himansu, P. C. E. Jorgenson, A modified space-time conservation element and solution element
method for Euler and Navier-Stokes equations. AIAA 99-3277 (1999) 10.

[[#citeF-38|[38]]] Y.-I. Lim, S.-C. Chang, S. B. Jørgensen, A novel partial differential algebraic equation (PDAE) solver: iterative space-time conservation
element/solution element (CE/SE) method. Computers & Chemical Engineering 28 (8) (2004) 1309–1324.

<div id=cite-39></div>
[[#citeF-39|[39]]] S.-T. J. Yu, L. Yang, R. L. Lowe, S. E. Bechtel, Numerical simulation of linear and nonlinear waves in hypoelastic solids by the CESE method. Wave Motion 47 (3) (2010) 168–182.

<div id=cite-40></div>
[[#citeF-40|[40]]] S.-C. Chang, X.-Y. Wang, C.-Y. Chow, The space-time conservation element and solution element method - a new high-resolution and
genuinely multidimensional paradigm for solving conservaton laws. Technical Memo TM 1998-208843, NASA, Lewis Research Center, NASA (August 1998).

<div id=cite-41></div>
[[#citeF-41|[41]]] M. Lukácová, K. W. Morton, G. Warnecke, Evolution Galerkin methods for hyperbolic systems in two space dimensions, Mathematics of
Computation 69 (232) (2000) 1355–1384.

<div id=cite-42></div>
[[#citeF-42|[42]]] Qurrat-Ul-Ain, S. Qamar, G. Warnecke, A high-resolution space-time conservative method for non-linear hyperbolic conservation laws.
International Journal of Computational Methods 4 (2) (2007) 223–247.

<div id=cite-43></div>
[[#citeF-43|[43]]] J. J. Stoker, Water Waves: The Mathematical Theory with Applications, Interscience Publishers, New York, 1957.

<div id=cite-44></div>
[[#citeF-44|[44]]] O. Delestre, C. Lucas, P.-A. Ksinant, F. Darboux, C. Laguerre, T.-N.-T. Vo, F. James, S. Cordier, SWASHES: a compilation of shallow water analytic solutions for hydraulic and environmental studies. International Journal for Numerical Methods in Fluids 72 (3) (2013) 269–300.

<div id=cite-45></div>
[[#citeF-45|[45]]] R. J. Fennema, M. H. Chaudhry, Implicit methods for two-dimensional unsteady free-surface flows. Journal of Hydraulic Research 27 (3)
(1989) 321–332.

<div id=cite-46></div>
[[#citeF-46|[46]]] C. Zoppou, S. Roberts, Numerical solution of the two-dimensional unsteady dam break. Applied Mathematical Modelling 24 (7) (2000) 457–475.

<div id=cite-47></div>
[[#citeF-47|[47]]] F. Alcrudo, P. Garcia-Navarro, A high-resolution Godunov-type scheme in finite volumes for 2D shallow-water equations. International
Journal for Numerical Methods in Fluids 16 (1993) 489–505.
</div>