% Template for submitting a manuscript to the Journal of Engineering of Catholic University of Petropolis. % version 2.0 dated 28 August 2012. % This file was adapted from cilamce 2011. Modifications may be freely made, % provided the edited file is saved under a different name %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% %% %% Important note on usage %% %% ----------------------- %% %% This file must be compiled with PDFLaTeX %% %% Using standard LaTeX will not work! %% %% %% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \documentclass[oneside,a4paper,portuguese,links]{EngUCP} % \usepackage{graphicx} \usepackage{amsmath,amsfonts} % \usepackage[portuguese,algoruled,boxed,linesnumbered,commentsnumbered]{algorithm2e} \usepackage{float} \title{COMPARAÇÃO DOS MODELOS DE DARCY, FORCHHEIMER E BARREE E CONWAY NA SIMULAÇÃO DE RESERVATÓRIOS DE GÁS} %\author[1]{First A. Author} %\author[1,2]{Second B. Author} %\author[3]{Third C. Author} % %% %\affil[1]{Department of Engineering and Computation, Catholic %University of Petrópolis , Petrópolis, 25685-070, RJ, Brazil} %% %\affil[2]{National Laboratory of Computing Science, Petrópolis, %25651-075, RJ, Brazil } % %\affil[3]{Department of Civil Engineering, Pontifical Catholic %University of Rio de Janeiro, Rio de Janeiro, 25651-075, RJ, %Brazil } %% NOTE: IF ALL AUTHORS BELONG TO THE SAME AFFILIATION %% USE THE `\voidaffil' MACRO FOR THE AFFILIATION CODE. %% Example: \author[\voidaffil]{Grazione de Souza} \author[\voidaffil]{Helio Pedro Amaral Souto} % \author[\voidaffil]{Third C. Author} % \author[\voidaffil]{Fourth D. Author} % % \affil[\voidaffil]{Instituto Politécnico, Universidade do Estado do Rio de Janeiro\\ Caixa Postal 97282, 28610-974 Nova Friburgo, RJ, Brasil} \pagestyle{empty} \begin{document} \vspace{3cm} \maketitle %E-mail Addresses: Put here! \let\oldthefootnote\thefootnote \renewcommand{\thefootnote}{\fnsymbol{footnote}} \footnotetext[1]{ E-mail addresses: \url{gsouza@iprj.uerj.br}, \url{helio@iprj.uerj.br} .} \let\thefootnote\oldthefootnote \begin{keywords} Equação de Forchheimer, Escoamento de Gás, Fluxo Não-Darcy, Método de Diferenças Finitas, Modelo de Barree e Conway. \end{keywords} \begin{abstract} Neste trabalho, um simulador numérico tridimensional de reservatórios de hidrocarbonetos foi utilizado com o intuito de comparar os resultados obtidos com diferentes modelos para a equação de {\it momentum}, e considerando o escoamento de gases em meios porosos heterogêneos. A conhecida lei de Darcy, válida para escoamentos a baixo número de Reynolds, foi aplicada e os efeitos inerciais foram modelados usando: i) a equação de Forchheimer e ii) o modelo de Barree e Conway. O conjunto de equações governantes, para o escoamento estudado, é formado pela conservação da massa, pelo balanço da quantidade de movimento linear e por uma equação de estado. Uma equação diferencial parcial não-linear para a variável pressão é obtida, a partir das equações de balanço e de estado, e é discretizada via o método de Diferenças Finitas com uma formulação implícita no tempo. O sistema de equações algébricas resultante do processo de dicretização é resolvido numericamente por um método do tipo Gradientes Conjugados com pré-condicionamento, apropriado para matrizes simétricas e não-simétricas. As soluções, obtidas para a pressão no poço produtor de gás, mostraram a influência das diferentes equações de {\it momentum} na modelagem do escoamento e reforçaram a necessidade de uma caracterização de reservatórios acurada, para que modelos apropriados para o uso em escoamentos inerciais sejam aplicados no gerenciamento da produção de hidrocarbonetos. \end{abstract} % -------------------------------------------------------------------------- \section{INTRODUÇÃO} % -------------------------------------------------------------------------- Desde a década de 1960, a simulação numérica de reservatórios de petróleo vem aumentado a sua importância na indústria de petróleo e gás. Simulações numéricas têm sido utilizadas no estudo, por exemplo, de escoamentos em meios porosos, da gestão de reservatórios, da dinâmica de poços e da caracterização de reservatórios. Portanto, a modelagem numérica vem sendo empregada, com êxito, a fim de se obter um processo otimizado de produção nos poços, de modo a se garantir que a recuperação de hidrocarbonetos seja maximizada. Conforme já mencionado, a simulação numérica de reservatórios, portadores de hidrocarbonetos, tem se tornado cada vez mais relevante nas atividades de gerenciamento de reservatórios. Para alguns casos mais simples, como o escoamento unidimensional de um fluido ligeiramente compressível, existem soluções analíticas para a determinação do fluxo em meios porosos. Por outro lado, para sistemas mais complexos e escoamentos mais realísticos, a simulação numérica tem se mostrado uma ferramenta valiosa no estudo do escoamento de fluidos compressíveis e incompressíveis em meios porosos heterogêneos~\cite{Ertekin2001}. Em se tratando do escoamento de gases, as equações governantes são não-lineares devido às propriedades do fluido serem dependentes da pressão reinante no reservatório~\cite{He2013,Jmili2011}. Além disso, não-linearidades podem aparecer em função dos efeitos inerciais em escoamentos do tipo não-Darcy~\cite{Tiab2004}. Exemplos de problemas relevantes de simulações numéricas envolvendo o escoamento de gases podem ser em encontrados em~\cite{DeSouza2013,He2013,Jmili2011,Tavares2006}. Em geral, escoamentos em meios porosos são modelados a partir da lei de Darcy. Entretanto, em alguns casos, os efeitos inerciais podem levar a resultados que não estão de acordo com essa lei. Situações onde estes efeitos podem ser relevantes ocorrem no escoamento de gases em meios fraturados e no escoamento próximo aos poços. A equação de Forchheimer é um exemplo bem conhecido de um modelo para escoamentos inerciais~\cite{Tavares2006}. Barree e Conway~\cite{Barree2004} também propuseram um outro modelo que incorpora os efeitos inerciais em escoamentos em meios porosos. Segundo Barree e Conway~\cite{Barree2004}, normalmente, o fluido deve acelerar ao passar por gargantas de poro estreitas e desacelerar quando entra em poros maiores. Portanto, as repetidas mudanças de velocidade conduzem a quedas de pressão associadas aos efeitos inerciais~\cite{Barree2004}. Estes efeitos são responsáveis pela redução da permeabilidade ``aparente'' associada ao escoamento do tipo não-Darcy. O principal objetivo deste trabalho é o estudo comparativo das previsões fornecidas, para a simulação numérica do escoamento monofásico, em um reservatório tridimensional de gás com um poço produtor vertical, pela lei de Darcy tradicional, pela equação de Forchheimer e pelo modelo de Barree e Conway. % ------------------------------------------------------------------------ \section{MODELAGEM FÍSICO-MATEMÁTICA} % ------------------------------------------------------------------------ Uma equação diferencial parcial, para a determinação da pressão do escoamento no meio poroso, é obtida usando-se as equações de balanço de massa, quantidade de movimento e de estado, considerando-se as seguintes hipóteses: \begin{itemize} \item meio poroso heterogêneo; \item permeabilidades constantes; \item compressibilidade da rocha pequena e constante; \item fluido Newtoniano; \item reservatório de gás seco; \item fluido de composição química constante; \item escoamento monofásico e isotérmico. \end{itemize} % -------------------------------------------------------------------------- \subsection{Equações Governantes para o Escoamento em Meio Poroso} % -------------------------------------------------------------------------- A equação da conservação de massa, no escoamento em meios porosos, pode ser escrita como \begin{equation} \frac{\partial}{\partial t}\left(\phi\rho\right)+\nabla\cdot\left(\rho\mathbf{v}\right)-\dfrac{q_m}{V_b}=0, \label{eq:engres4.1} \end{equation} \noindent onde $\phi$ é a porosidade, $\rho$ a massa específica, $\mathbf{v}$ a velocidade superficial de escoamento, $q_m$ indica um termo fonte e $V_b$ é o volume da formação (rocha e poros). Para o fluxo de um gás real, a massa específica é uma função não-linear da pressão e introduz-se o Fator-Volume-Formação, simbolizado por $B$~\cite{Ertekin2001}, de modo que \begin{equation} B=\frac{V}{V_{sc}}=\frac{\rho_{sc}}{\rho}=\frac{p_{sc}TZ}{T_{sc}p}, \label{eq:ertekin2.83} \end{equation} \noindent onde $V$ é o volume de fluido nas condições de reservatório, $V_{sc}$ o volume de fluido em condições padrão, $T$ a temperatura, $p$ a pressão, $Z$ simboliza o fator de compressibilidade (função de $p$ e $T$), e $\rho_{sc}$, $p_{sc}$ e $T_{sc}$ são a massa específica, a pressão e a temperatura em condições padrão, respectivamente. Em seguida, uma forma geral, que expressa a conservação da quantidade de movimento, é escrita modificando-se a lei clássica de Darcy \begin{equation} \mathbf{v}=-\frac{\mathbf{k}_{mod}}{\mu}[\nabla p-\rho g\nabla D], \label{eq:darcy_kapp} \end{equation} \noindent onde $\mathbf{k}_{mod}=\mathbf{k}$ para a lei original de Darcy, $\mathbf{k}_{mod}=\boldsymbol{\delta}\mathbf{k}$ para a equação de Forchheimer e $\mathbf{k}_{mod}=\mathbf{k}_{BC}\mathbf{k}$ para o modelo de Barree e Conway. O tensor $\mathbf{k}$ representa o tensor de permeabilidade absoluta do meio poroso, $\mu$ a viscosidade do fluido, $\nabla D$ o gradiente de profundidade, $g$ a magnitude da aceleração da gravidade e $\boldsymbol{\delta}$ e $\mathbf{k}_{BC}$ serão introduzidos na sequência, em função do modelo inercial. Para o escoamento de um gás, a viscosidade deve ser modelada como sendo uma função da temperatura e da pressão~\cite{Ertekin2001}. Após a realização de experimentos com escoamentos de gases, Forchheimer modificou a lei de Darcy com o objetivo de incluir os efeitos inerciais no balanço da quantidade de movimento~\cite{Zeng2008} \begin{equation} \boldsymbol{\delta}=\left(\mathbf{1}+\frac{\rho\boldsymbol{\beta} \mathbf{k} |\mathbf{v}|}{\mu}\right)^{-1}, \label{eq:deviation} \end{equation} \noindent na qual $\boldsymbol{\beta}$ é descrito na forma \begin{equation} \boldsymbol{\beta}={\displaystyle c_{1}\mathbf{k}^{c_{2}}\phi^{c_{3}}\boldsymbol{\tau}_{p}^{c_{4}},} \label{eq:eq5.27-1} \end{equation} \noindent onde $\boldsymbol{\tau}_{p}$ representa o tensor de tortuosidade do meio poroso, enquanto que $c_{1}$, $c_{2}$, $c_{3}$ e $c_{4}$ são constantes para uma dada formação rochosa. $\boldsymbol{\beta}$ pode ser obtido através de dados de produção ou a partir de experimentos de fluxo em meios porosos~\cite{Li2001}. Barree e Conway~\cite{Barree2004} propuseram uma outra forma de incorporar os efeitos inerciais na modelagem do escoamento em meios porosos, na qual a permeabilidade é considerada de forma diferenciada~\cite{Lai2009}. A permeabilidade $\mathbf{k}$ é substituída por uma permeabilidade aparente $\mathbf{k}_{app}=\mathbf{k}_{BC}\mathbf{k}$, dependente da vazão e com $\mathbf{k}_{BC}$ sendo definido por \begin{equation} \mathbf{k}_{BC}=\left[\mathbf{k}_{mr}+\dfrac{\mathbf{1}-\mathbf{k}_{mr}}{1+\left(\displaystyle\frac{\rho |\mathbf{v}|}{\mu\boldsymbol{\tau}_c}\right)}\right] , \label{eq:k_app} \end{equation} \noindent onde $\mathbf{k}_{mr}=\mathbf{k}_{min}\mathbf{k}^{-1}$, sendo $\mathbf{k}_{min}$ um valor mínimo da permeabilidade para altas vazões de escoamento e $\boldsymbol{\tau}_{c}$ possui a dimensão do inverso de um comprimento característico~\cite{Barree2004}. De acordo com Lai~et al.~\cite{Lai2009}, o modelo de Barree e Conway apresenta algumas vantagens quando comparado a outras equações aplicadas na modelagem do balanço de momento linear no escoamento em meios porosos. Usando o modelo de Barree e Conway, uma única equação pode ser usada para descrever uma faixa inteira de escoamentos. Além disso, regiões de transição entre os regimes de escoamento seriam melhor capturadas. Lai~et al.~\cite{Lai2009} também relataram que o modelo de Barree e Conway provê um {\it plateau} relacionado às altas vazões de escoamento (descrito na literatura) e também recupera resultados da lei clássica de Darcy e da equação de Forchheimer, para algumas faixas de velocidades de escoamento. Considerando-se que os modelos de Forchheimer e Barree e Conway levassem aos mesmos resultados, isto resultaria em \begin{equation} \left(\mathbf{1}+\frac{\rho\boldsymbol{\beta} \mathbf{k} |\mathbf{v}|}{\mu}\right)^{-1}=\mathbf{k}_{mr}+\dfrac{\mathbf{1}-\mathbf{k}_{mr}}{1+\left(\displaystyle\frac{\rho |\mathbf{v}|}{\mu\boldsymbol{\tau}_c}\right)}. \label{eq:forchheimer_equal_barree} \end{equation} Jones~\cite{Jones1987} baseado em medidas obtidas para $\boldsymbol{\beta}$, para diferentes valores de permeabilidade, sugeriu uma relação entre $\boldsymbol{\beta}$ e $\mathbf{k}$ na forma \begin{equation} \boldsymbol{\beta}={\displaystyle 6,15 \times 10^{10}\,\mathbf{k}^{-1,55},} \label{eq:Jones1987} \end{equation} \noindent para $\mathbf{k}$ expresso em md e valores de $\boldsymbol{\beta}$ dados em 1/ft, e ela é conhecida na literatura como a correlação de Jones. De fato, existem diversas propostas na literatura para a determinação de $\boldsymbol{\beta}$, baseadas em cálculos teóricos, em experimentos em laboratório ou em dados de produção. Uma revisão com diferentes correlações para a determinação de $\boldsymbol{\beta}$, indicadas para o seu uso em diferentes aplicações, pode ser encontrada em Li e Engler~\cite{Li2001a}. Dentre as correlações revisadas pelos autores encontram-se \begin{equation} \boldsymbol{\beta}=1,82\times10^8\, \mathbf{k}^{-5/4}\phi^{-3/4} \label{eq:Janicek1955} \end{equation} \noindent proposta por Janicek e Katz~\cite{Janicek1955}, na qual $\mathbf{k}$ é dado em md e $\boldsymbol{\beta}$ em 1/cm, e a correlação \begin{equation} \boldsymbol{\beta}=8,91\times10^8\, \mathbf{k}^{-1}\phi^{-1}\boldsymbol{\tau_p} \label{eq:Liu1995} \end{equation} \noindent sugerida por Liu et al.~\cite{Liu1995}, com $\mathbf{k}$ em md e $\mathbf{\beta}$ em 1/ft. Os valores dos coeficientes usados na Eq.~\eqref{eq:eq5.27-1}, por cada um dos três modelos apresentados, encontram-se agrupados na Tabela~\ref{table1}. \begin{table}[t] \centering \caption{Diferentes correlações para o cálculo de $\beta$.}\label{table1}\vspace{0.25cm} \begin{tabular}{ccccc} \hline Modelo &$c_1$ &$c_2$ &$c_3$ &$c_4$ \\ \hline Jones~\cite{Jones1987} &$6,15 \times 10^{10}$ &$-1,55$ &$0$ &$0$ \\ Janicek e Katz~\cite{Janicek1955} &$1,82 \times 10^{8}$ &$-5/4$ &$-3/4$ &$0$ \\ Liu et al.~\cite{Liu1995} &$8,91 \times 10^{8}$ &$-1$ &$-1$ &$1$ \\ \hline \end{tabular} \end{table} Barree e Conway~\cite{Barree2004} relataram que os resultados obtidos com a correlação proposta por Jones~\cite{Jones1987} desviam dos valores experimentais para faixas de permeabilidade elevadas, e que neste caso, uma função não-linear se ajustaria melhor aos dados. Considerando as faixas de $\boldsymbol{\beta}$ estudadas em~\cite{Jones1987}, para as quais os desvios com relação à equação de Forchheimer são pequenos, $\boldsymbol{\tau}_c$ pode ser determinado usando-se \begin{equation} \boldsymbol{\tau}_c={\displaystyle \mathbf{k}^{-1}\boldsymbol{\beta}^{-1},} \label{eq:tauc_beta} \end{equation} \noindent o qual é obtido tomando-se $\mathbf{k}_{mr}=\mathbf{0}$ na Eq.~\eqref{eq:forchheimer_equal_barree}. % -------------------------------------------------------------------------- \subsection{Equação da Difusividade Hidráulica} % -------------------------------------------------------------------------- Com o objetivo de obter-se uma equação para a variável primitiva pressão, o balanço de massa, a equação de estado para um gás real e a expressão modificada para o balanço da quantidade de movimento linear, Eq.~\eqref{eq:darcy_kapp}, são consideradas como equações fundamentais. Primeiramente, substitui-se a Eq.~\eqref{eq:darcy_kapp} na Eq.~\eqref{eq:engres4.1}, aplicando-se em conjunto a Eq.~\eqref{eq:ertekin2.83}. O resultado é, então, multiplicado por $V_b$. Neste desenvolvimento considera-se que \begin{equation} q_{m}=-J_w(p-p_{wf})\rho_{sc}, \end{equation} \noindent onde $J_w$ é o índice de produtividade e $p_{wf}$ a pressão no poço produtor~\cite{Ertekin2001} e utiliza-se a relação \begin{equation} \phi=\phi^{0}\left[1+c_{\phi}\left(p-p^{0}\right)\right], \end{equation} \noindent na qual $\phi^{0}$ e $p^{0}$ são a porosidade e a pressão de referência e $c_{\phi}$ a compressibilidade da rocha~\cite{Ertekin2001}. Após realizadas todas as devidas substituições é possível escrever-se \begin{equation} \Gamma\dfrac{\partial p}{\partial t}-V_b\nabla\cdot\left[\dfrac{\mathbf{k}_{mod}}{B\mu}\nabla p\right]+J_w(p-p_{wf})+\Lambda_G=0, \label{eq:engres4.1_5} \end{equation} \noindent onde \begin{equation} \Gamma=V_b\left[\phi^0\left(\frac{c_{\phi}}{B}\right)+\phi \frac{\partial}{\partial p}\left(\frac{1}{B}\right)\right], \end{equation} \noindent sendo $\Lambda_G=V_bG$, com $G$ incluindo os efeitos gravitacionais. O índice $J_w$ é uma função de alguns parâmetros, tais como das propriedades do fluido, da rocha e da geometria do reservatório~\cite{Dumkwu2012}. Considerando-se o sistema de coordenadas Cartesianas, a condição inicial para a pressão é dada por \begin{equation} p(x,y,z,t=0)=p_{ini}(z), \label{eq:EDH_delta3_condicao_inicial} \end{equation} \noindent onde $p_{ini}$ é uma função de $z$, a profundidade do reservatório. Em geral, é dada a pressão $P_i$ no topo do reservatório e a distribuição $p_{ini}$ é determinada pela hidrostática. As condições de fronteira externa para o reservatório são do tipo fluxo nulo \begin{equation} \left(\dfrac{\partial p}{\partial x}\right)_{x=0,L_{x}}=\left(\dfrac{\partial p}{\partial y}\right)_{y=0,L_{y}}=\left(\dfrac{\partial p}{\partial z}\right)_{z=0,L_{z}}=0, \label{eq:EDH_cilindrica_condicao_contorno_externa} \end{equation} \noindent onde $L_x$, $L_y$ e $L_z$ são os comprimentos do reservatório nas direções $x$, $y$ e $z$, respectivamente. % -------------------------------------------------------------------------- \section{SOLUÇÃO NUMÉRICA} % -------------------------------------------------------------------------- Na discretização da equação diferencial parcial~\eqref{eq:engres4.1_5}, escrita em termos da pressão, é empregado o método de Diferenças Finitas~\cite{Ertekin2001}, enquanto que o método dos Gradientes Conjugados pré-condicionado (PCG) é utilizado na resolução numérica do sistema de equações algébricas. % -------------------------------------------------------------------------- \subsection{Discretização da Equação de Difusividade Hidráulica} % -------------------------------------------------------------------------- A Fig.~\ref{fig:Reservatorio3d} mostra uma representação esquemática da discretização do reservatório, usando o sistema tridimensional de coordenadas Cartesianas. Os índices $i$, $j$ e $k$ indicam os centros das células ao longo das direções $x$, $y$ e $z$, respectivamente. Os índices fracionários $i\pm 1/2,j,k$; $i,j\pm 1/2,k$; e $i,j,k\pm 1/2$ representam as fronteiras das células. \begin{figure}[h!] \centering \includegraphics[scale=0.47]{disc2.png} \caption{Discretização do reservatório empregando uma malha estruturada.} \label{fig:Reservatorio3d} \end{figure} Neste trabalho, supõe-se um meio anisotrópico com $\mathbf{k}_{mod}$ sendo um tensor diagonal, com componentes não nulas $k_{mod,xx}$, $k_{mod,yy}$, e $k_{mod,zz}$. Portanto, \begin{equation} \Gamma\dfrac{\partial p}{\partial t} -\dfrac{\partial}{\partial x}\!\left(\dfrac{k_{mod,xx}A_x}{B\mu}\dfrac{\partial p}{\partial x}\right)\!\!\Delta x -\dfrac{\partial}{\partial y}\!\left(\dfrac{k_{mod,yy}A_y}{B\mu}\dfrac{\partial p}{\partial y}\right)\!\!\Delta y -\dfrac{\partial}{\partial z}\!\left(\dfrac{k_{mod,zz}A_z}{B\mu}\dfrac{\partial p}{\partial z}\right)\!\!\Delta z -q_{sc}\!=\!-\Lambda_G, \label{eq:engres4.1_6} \end{equation} \noindent considerando-se que $q_{sc}=q_{m}/\rho_{sc}$, $V_b=\Delta x \Delta y \Delta z$, $A_x=\Delta y \Delta z$, $A_y=\Delta x \Delta y$ e $A_z=\Delta x \Delta y$. Usando-se blocos centrados, aproximações do tipo diferença avançada no tempo e centradas no espaço~\cite{Ertekin2001} chega-se à forma final discretizada para a formulação implícita \[ \dfrac{\Gamma^{n+1}_{i,j,k}}{\Delta t}\left(p_{i,j,k}^{n+1}-p_{i,j,k}^{n}\right) -T^{n+1}_{x,i-\frac{1}{2},j,k}\left(p_{i-1,j,k}^{n+1}-p_{i,j,k}^{n+1}\right) -T^{n+1}_{x,i+\frac{1}{2},j,k}\left(p_{i+1,j,k}^{n+1}-p_{i,j,k}^{n+1}\right) \] \[ -T^{n+1}_{y,i,j-\frac{1}{2},k}\left(p_{i,j-1,k}^{n+1}-p_{i,j,k}^{n+1}\right) -T^{n+1}_{y,i,j+\frac{1}{2},k}\left(p_{i,j+1,k}^{n+1}-p_{i,j,k}^{n+1}\right) \] \[ -T^{n+1}_{z,i,j,k-\frac{1}{2}}\left(p_{i,j,k-1}^{n+1}-p_{i,j,k}^{n+1}\right) -T^{n+1}_{z,i,j,k+\frac{1}{2}}\left(p_{i,j,k+1}^{n+1}-p_{i,j,k}^{n+1}\right) \] \begin{equation} -J^{n+1}_{w,i,j,k}\left(p_{wf,i,j,k}^{n+1}-p_{i,j,k}^{n+1}\right)=-\left(\Lambda_G\right)^{n+1}_{i,j,k} \label{discrete_xyz} \end{equation} \noindent onde as transmissibilidades foram introduzidas como sendo dadas por \begin{equation} T^{n+1}_{x,i\pm\frac{1}{2},j,k}\!=\!\left(\dfrac{k_{mod,xx}A_x}{B\mu\Delta x}\right)^{n+1}_{i\pm\frac{1}{2},j,k}; \ \label{eq:transmissibility_x} \end{equation} \begin{equation} T^{n+1}_{y,i,j\pm\frac{1}{2},k}\!=\!\left(\dfrac{k_{mod,yy}A_y}{B\mu\Delta y}\right)^{n+1}_{i,j\pm\frac{1}{2},k}; \ \label{eq:transmissibility_y} \end{equation} \begin{equation} T^{n+1}_{z,i,j,k\pm\frac{1}{2}}\!=\!\left(\dfrac{k_{mod,zz}A_z}{B\mu\Delta z}\right)^{n+1}_{i,j,k\pm\frac{1}{2}} \label{eq:transmissibility_z} \end{equation} \noindent sendo que uma média harmônica é utilizada na determinação das transmissibilidades nas interfaces das células e~\cite{Ertekin2001} \begin{equation} \Gamma^{n+1}_{i,j,k}=V_{b,i,j,k}\left\{\dfrac{\phi^{0}c_{\phi}}{B^{n+1}}+\dfrac{\phi^{n}}{B^{n}}\left[\dfrac{\left(B^{n}/B^{n+1}\right)-1}{p^{n+1}-p^{n}}\right]\right\} _{i,j,k}. \label{gamma_xyz} \end{equation} % -------------------------------------------------------------------------- \subsection{Resolução Numérica em Termos da Pressão} % -------------------------------------------------------------------------- Escrita para cada célula, a Eq.~\eqref{discrete_xyz} leva a um conjunto de equações algébricas não-lineares em termos da variável pressão. A fim de que se possa empregar técnicas de resolução de sistemas de equações lineares, necessita-se linearizar este sistema. Admitindo-se que para uma iteração da transmissibilidade $T^{n+1}\cong T^{n+1,v}$~\cite{Ertekin2001}, a Eq.~\eqref{discrete_xyz} pode ser reescrita como \[ \left(\frac{\Gamma^{n+1,v}_{i,j,k}}{\Delta t}\right)\left(p_{i,j,k}^{n+1,v+1}-p_{i,j,k}^{n}\right) -T^{n+1,v}_{x,i-\frac{1}{2},j,k}\left(p_{i-1,j,k}^{n+1,v+1}-p_{i,j,k}^{n+1,v+1}\right) -T^{n+1,v}_{x,i+\frac{1}{2},j,k}\left(p_{i+1,j,k}^{n+1,v+1}-p_{i,j,k}^{n+1,v+1}\right) \] \[ -T^{n+1,v}_{y,i,j-\frac{1}{2},k}\left(p_{i,j-1,k}^{n+1,v+1}-p_{i,j,k}^{n+1,v+1}\right) -T^{n+1,v}_{y,i,j+\frac{1}{2},k}\left(p_{i,j+1,k}^{n+1,v+1}-p_{i,j,k}^{n+1,v+1}\right) \] \[ -T^{n+1,v}_{z,i,j,k-\frac{1}{2}}\left(p_{i,j,k-1}^{n+1,v+1}-p_{i,j,k}^{n+1,v+1}\right) -T^{n+1,v}_{z,i,j,k+\frac{1}{2}}\left(p_{i,j,k+1}^{n+1,v+1}-p_{i,j,k}^{n+1,v+1}\right) \] \begin{equation} -J^{n+1,v}_{w,i,j,k}\left(p_{wf,i,j,k}^{n+1,v+1}-p_{i,j,k}^{n+1,v+1}\right)=-\left(\Lambda_G\right)^{n+1,v}_{i,j,k}, \label{discrete_xyz1} \end{equation} \noindent sendo os níveis de iteração indicados por $v$ (valores conhecidos) e $v+1$ (valores desconhecidos). O mesmo procedimento é usado para os termos $\Gamma_{i,j,k}$ e $J_w$. Para descrever-se o acoplamento poço-reservatório para um escoamento monofásico, neces\-si\-ta-se de uma equação para a determinação da taxa de escoamento no poço e precisa-se calcular $J_w$. Neste trabalho, somente poços verticais são tratados e os efeitos de perdas devidas ao atrito e os efeitos inerciais são desconsiderados no interior do poço. Na modelagem do acoplamento, a vazão total de fluxo no poço, $Q_{sc}$, deve ser igual ao somatório do fluxo de gás deixando as células (em direção ao poço) que estão conectadas ao poço, assim~\cite{Ertekin2001} \begin{equation} Q_{sc}=-\sum^{n_z}_{k=1}J^{n+1,v}_{w,i,j,k}\left[p^{n+1,v+1}_{i,j,k}-\left(p_{wf}\right)^{n+1,v+1}_{i,j,k}\right] \label{Jw3} \end{equation} \noindent para um poço trespassando $n_z$ camadas. Adotando-se um nível como sendo o de referência ($k=W$) para o cômputo da pressão no poço e levando-se em consideração os efeitos gravitacionais no interior do mesmo, é possível escrever~\cite{Ertekin2001} \begin{equation} \left(p_{wf}\right)_{i,j,W}^{n+1,v+1}=\dfrac{\displaystyle\sum^{n_z}_{k=1}J^{n+1,v}_{w,i,j,k}\left[ \left(p_{wf}\right)^{n+1,v+1}_{i,j,k}-\left(\rho g\right)^{n+1,v}_{i,j,k+1/2}\left(z_{i,j,k}-z_{i,j,W}\right)\right]+Q_{sc}}{\displaystyle\sum^{n_z}_{k=1}J^{n+1,v}_{w,i,j,k}} \label{Jw4} \end{equation} \noindent para $\rho^{n+1,v}_{i,j,k+1/2}$ sendo calculado como uma média entre os valores $\rho^{n+1,v}_{i,j,k}$ e $\rho^{n+1,v}_{i,j,W}$, e \begin{equation} \left(p_{wf}\right)_{i,j,k}^{n+1,v+1}=\left(p_{wf}\right)_{i,j,W}^{n+1,v+1}+\left(\rho g\right)^{n+1,v}_{i,j,k+1/2}\left(z_{i,j,k}-z_{i,j,W}\right). \label{Jw5} \end{equation} Para o índice de produtividade, $J_w$, usa-se \begin{equation} J_{w,i,j,k}^{n+1,v}=\left[\frac{2\pi\sqrt{k_{mod,xx}k_{mod,yy}}\Delta z}{\mu B\ln\left(r_{eq}/r_{w}\right)}\right]_{i,j,k}^{n+1,v} \label{Jw2} \end{equation} \noindent e \begin{equation} r_{eq,i,j,k}^{n+1,v}=0,28\left[\frac{\sqrt{k_{mod,yy}/k_{mod,xx}}\left(\Delta x\right)^{2}+\sqrt{k_{mod,xx}/k_{mod,yy}}\left(\Delta y\right)^{2}}{\left(k_{mod,yy}/k_{mod,xx}\right)^{1/4}+\left(k_{mod,xx}/k_{mod,yy}\right)^{1/4}}\right]_{i,j,k}^{n+1,v}, \label{req} \end{equation} \noindent onde $r_{eq}$ é o raio equivalente~\cite{Peaceman1983} e $r_w$ o raio do poço. A fim de que os efeitos inerciais sejam levados em consideração, as permeabilidades modificadas são utilizadas na determinação de $J_w$ e $r_{eq}$ empregando-se as Eqs.~\eqref{Jw2} e~\eqref{req}, respectivamente~\cite{Al_Mohannadi2004}. As Eqs.~\eqref{discrete_xyz1} e~\eqref{Jw4} formam um sistema de equações para a determinação da pressão no reservatório e no poço. Este sistema de equações, linearizado, é resolvido numericamente por uma técnica de fatorização aproximada. Esta técnica possui a capacidade de tratar de problemas que representam situações práticas de interesse, nas quais as heterogeneidades estão presentes e, devido às escalas consideradas, podem impactar na variabilidade das propriedades da matriz dos coeficientes (condição da matriz). A versão para matrizes simétricas e assimétricas do método dos Gradientes Conjugados, com pré-condicionamento, foi empregada para a resolução do sistema de equações algébricas~\eqref{discrete_xyz1} e~\eqref{Jw4}, utilizando uma matriz de precondicionamento de Jacobi (diagonal)~\cite{Ertekin2001,Saad2003}. A solução do processo iterativo é alcançada quando \begin{equation} max\left|p^{n+1,v+1}-p^{n+1,v}\right|10^{-1}$ dia). Em função dos resultados apresentados após as primeiras horas de produção de gás, os valores $n_x=n_y=81$ foram escolhidos para serem utilizados em todas as demais simulações. Os resultados do simulador numérico foram validados mediante a comparação deles com os respectivos valores produzidos por~\cite{DeSouza2013}, para as mesmas condições de escoamento e a lei clássica de Darcy, e sem considerar os efeitos do fluxo dentro do poço de produção. As curvas utilizadas na comparação encontram-se na~Fig.~\ref{fig:validation} e pode-se constatar que os valores deste trabalho correspondem muito bem aos obtidos com o simulador numérico de~\cite{DeSouza2013}, cuja funcionalidade e acurácia foram testadas via uma confrontação direta com o simulador comercial~IMEX~\cite{CMG2009}. Deve-se destacar que uma muito boa concordância foi verificada entre os resultados produzidos pelos código numérico IMEX e o simulador de de~Souza~\cite{DeSouza2013}. \begin{figure}[h] \centering \includegraphics[scale=0.5]{validation2.png} \caption{Validação do simulador numérico.} \label{fig:validation} \end{figure} \subsection{Os Modelos Inerciais} Utilizando-se os valores previamente introduzidos na Tabela~\ref{table2}, simulações numéricas foram realizadas e comparados os resultados oriundos da aplicação dos diferentes modelos inerciais. No caso da equação de Forchheimer, considerando a correlação de Jones~\cite{Jones1987} (Eq.~\eqref{eq:Jones1987}), os valores de $c_1$, $c_2$ e $c_3$ e $c_4$ utilizados encontram-se na Tabela~\ref{table1} enquanto que para o modelo de Barree e Conway empregou-se $k_{mr}=0,5$ e $\tau_c=3,281$ m$^{-1}$. As Figs.~\ref{fig:q25}-\ref{fig:q100} mostram os resultados para a pressão no poço usando-se as três equações para o balanço da quantidade de movimento linear: i) lei clássica de Darcy ($\mathbf{k}_{mod}=\mathbf{k}$), ii) equação de Forchheimer ($\mathbf{k}_{mod}=\boldsymbol{\delta}\mathbf{k}$) e iii) modelo de Barree e Conway ($\mathbf{k}_{mod}=\mathbf{k}_{BC}\mathbf{k}$); tendo sido variada a vazão de produção $Q_{sc}$. Conforme o esperado, a equação de Forchheimer e o modelo de Barree e Conway conduziram a maiores quedas da pressão do poço, em função da captura dos efeitos de inércia. Também é possível observar que os efeitos tornam-se mais significativos à medida que aumenta-se a vazão, o que é fisicamente correto. Entretanto, na sequência das figuras, é interessante observar-se que a taxa de aumento da queda de pressão não é a mesma para os dois modelos (ao menos para o conjunto de parâmetros utilizados nas simulações), Forchheimer e Barree e Conway, que levam em consideração os efeitos inerciais e ambos diferem apreciavelmente da taxa de queda de pressão prevista pela lei original de Darcy. Tal constatação fica claramente evidenciada na Fig.~\ref{fig:q100}, quando a queda de pressão associada à equação de Forchheimer é maior do que àquela computada pelo modelo de Barree e Conway. Nas Figs.~\ref{fig:q25} e~\ref{fig:q50}, verifica-se que as variações de pressão, com relação à pressão inicial, foram mais significativas para o modelo de Barree e Conway. \begin{figure}[H] \centering \includegraphics[scale=0.5]{q25.png} \caption{Resultados para os diferentes modelos para equação de {\it momentum} e $Q_{sc}=$ -7,1591 $\times 10^4$ std m$^3$/dia.} \label{fig:q25} \end{figure} \begin{figure}[H] \centering \includegraphics[scale=0.5]{q50.png} \caption{Resultados para os diferentes modelos para equação de {\it momentum} e $Q_{sc}=$ -14,3182 $\times 10^4$ std m$^3$/dia.} \label{fig:q50} \end{figure} \begin{figure}[H] \centering \includegraphics[scale=0.5]{q100.png} \caption{Resultados para os diferentes modelos para equação de {\it momentum} e $Q_{sc}=$ -28,6364 $\times 10^4$ std m$^3$/dia.} \label{fig:q100} \end{figure} Em um outro teste, diferentes valores de $\tau_{c}$ foram utilizados para o modelo de Barree e Conway. Os resultados foram comparados com aqueles já obtidos com a lei original de Darcy e a equação de Forchheimer, para o poço produzindo a uma vazão constante de $Q_{sc}=$ $-28,6364\times 10^4$ std m$^3$/dia. O objetivo destas simulações foi o de quantificar a variação da perda de pressão para o modelo de Barree e Conway, para uma mesma vazão, ao se alterar um de seus parâmetros. Observando-se a Fig.~\ref{fig:tau}, percebe-se que o modelo de Barree e Conway pode conduzir a variações de pressão superiores ou inferiores às associadas à equação de Forchheimer, em função da variação do parâmetro $\tau_c$. Para valores mais baixos de $\tau_c$, as quedas de pressão foram mais elevadas, sendo esse resultado fisicamente compatível com a modelagem do fenômeno, uma vez que os efeitos inerciais são mais significativos neste caso. Finalmente, os resultados obtidos com o uso da equação de Forchheimer e diferentes correlações para a determinação de $\beta$ foram comparados: i)~Jones~\cite{Jones1987}, ii)~Janicek e Katz~\cite{Janicek1955} e iii)~Liu et al.~\cite{Liu1995}, com $\tau_p=\sqrt{2}$~\cite{Civan1998}. Ressalta-se, novamente, que os valores dos coeficientes $c_1$, $c_2$, $c_3$ e $c_4$ (da Eq.~\eqref{eq:eq5.27-1}) adotados pelos três modelos encontram-se disponíveis na Tabela~\ref{table1}. Uma confrontação dos valores de pressão, calculados para estes $\beta$, também é feita com aqueles já obtidos com a lei de Darcy (sem considerar os efeitos inerciais) e com o modelo de Barree e Conway, para $\tau_c$= 3,281 m$^{-1}$. Nas simulações, o poço produziu a uma vazão constante de $Q_{sc}=$ -28,6364 $\times 10^4$ std m$^3$/dia. Novamente, buscou-se quantificar a variação da queda de pressão, agora para a equação de Forchheimer, para uma mesma vazão e alterando um dos parâmetros utilizados na modelagem do escoamento. As curvas de pressão no poço para este conjunto de testes podem ser vistas na Fig.~\ref{fig:forchs2}. Nota-se que para uma mesma permeabilidade, diferentes correlações conduzem a resultados com diferentes desvios em relação à lei clássica de Darcy, podendo, em alguns casos, obter-se valores que encontram-se próximos dos resultados determinados com o modelo de Barree e Conway. Em função da grande quantidade de parâmetros que devem ser ajustados para cada um dos modelos inerciais, conforme visto, variações significativas dos valores da pressão podem ser encontradas. Portanto, a aplicabilidade dos mesmos depende de uma validação com resultados determinados experimentalmente ou com valores obtidos através da análise de dados de produção. \begin{figure}[H] \centering \includegraphics[scale=0.5]{tau.png} \caption{Resultados para os diferentes modelos para equação de {\it momentum}, incluindo a variação de $\tau_c$ no modelo de Barree e Conway (em todos os resultados, $Q_{sc}=$ -28,6364 $\times 10^4$ std m$^3$/dia).} \label{fig:tau} \end{figure} \begin{figure}[h] \centering \includegraphics[scale=0.5]{forchs2.png} \caption{Resultados para os diferentes modelos para equação de {\it momentum}, incluindo a variação da correlação para o termo $\beta$ utilizado na equação de Forchheimer (em todos os resultados, $Q_{sc}=$ -28,6364 $\times 10^4$ std m$^3$/dia).} \label{fig:forchs2} \end{figure} % -------------------------------------------------------------------------- \section{CONCLUSÕES}\label{Conclusao} % -------------------------------------------------------------------------- Neste artigo, empregou-se os modelos inerciais de Forchheimer e Barree e Conway no estudo da simulação numérica do escoamento de gás em reservatórios de hidrocarbonetos. O método de Diferenças Finitas e uma formulação implícita foram utilizados a fim de discretizar-se a equação diferencial parcial, escrita em termos da pressão e oriunda da equação de conservação de massa, da lei de Darcy modificada e da equação de estado para o gás real. Deste modo, pôde-se obter os valores numéricos para a pressão no poço ($p_{wf}$) considerando-se uma vazão ($Q_{sc}$) constante de produção. A partir da comparação dos modelos inerciais foi possível evidenciar diferenças significativas entre eles e capturar a influência dos efeitos inerciais na taxa de variação da queda de pressão no poço produtor. Em função dos resultados deste trabalho, evidenciou-se a necessidade de uma análise cuidadosa quando da simulação de reservatórios, devido à influência de parâmetros, tais como $\beta$ e $\tau_c$, nas previsões dos modelos de Forchheimer e de Barree e Conway. Na presença de escoamentos do tipo não-Darcy, os efeitos inerciais podem aumentar de modo apreciável a queda de pressão no poço produtor e é possível que, sob certas condições, diferentes modelos que incorporam os efeitos inerciais conduzam a resultados próximos uns dos outros. Isto reforça a necessidade de uma acurada caracterização de reservatórios, para que os resultados das previsões de produção sejam os mais fiéis à realidade, a partir da escolha adequada da equação de {\it momentum}. % -------------------------------------------------------------------------- \subsection*{\textit{Acknowledgments}} % -------------------------------------------------------------------------- \noindent Este trabalho contou com o apoio do Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq), através do processo 305958/2012-7, e da Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (CAPES). % \bibliographystyle{plainnat} %\bibliography{bibliography/ref} \bibliography{reucp_grazione} \end{document}