Ajuste de parâmetros de transporte de solutos no solo utilizando Matlab 6.5

Adjustment of soil solute transport parameters with Matlab 6.5

ARTIGOS CIENTÍFICOS

ENGENHARIA DE ÁGUA E SOLO

Ajuste de parâmetros de transporte de solutos no solo utilizando Matlab 6.5

Adjustment of soil solute transport parameters with Matlab 6.5

Anderson L. de Souza^I; Edson E. Matsura^II; Jarbas H. de Miranda^III; Alberto Colombo^IV

^IProfessor Assistente II da faculdade de tecnologia de Tatuí-SP

^IIProf. Titular da Faculdade de Engenharia Agrícola- Universidade de Campinas

^IIIProf. Associado do Depto. de Engenharia de Biossistemas da Escola Superior de Agricultura "Luiz de Queiroz" - ESALQ/USP

^IVProf. Associado, PhD, Departamento de Engenharia, UFLA. [e-mail]

RESUMO

O sucesso na utilização de modelos matemáticos no estudo do transporte de íons no solo está intimamente ligado à precisão com que os parâmetros de transporte envolvidos neste processo são estabelecidos. De maneira geral, tais parâmetros são determinados mediante a resolução de um problema de otimização não linear em que os dados experimentais, obtidos em ensaios de deslocamento miscível, são ajustados a um modelo teórico. Neste sentido, a utilização de softwares de alta performance no ajuste destes parâmetros mostra-se vantajosa, uma vez que, além da consistência e da disponibilidade de ferramentas numéricas preexistentes, possibilita a incorporação de novas rotinas de acordo com o fenômeno que se queira simular. Sendo assim, este trabalho teve como objetivo desenvolver, em ambiente MATLAB 6.5, uma rotina computacional para a otimização dos seguintes parâmetros de transporte: fator de retardamento (R) e coeficiente de dispersão (D). A rotina desenvolvida foi aplicada a dados experimentais de três ensaios de deslocamento miscível do íon potássio em colunas preenchidas com um Latossolo Vermelho- -Amarelo, fase arenosa. A qualidade dos ajustes obtidos foi avaliada utilizando-se do coeficiente de exatidão. Concluiu-se que a rotina proposta apresentou ótimo desempenho, o que, além de reforçar a consistência do método numérico utilizado, indica que a rotina proposta neste trabalho pode contribuir com o avanço dos estudos teóricos da dinâmica da água e de solutos em meios porosos não saturados.

Palavras-chave: curva de distribuição de efluentes, rotina computacional, dinâmica de solutos.

ABSTRACT

The successful use of mathematical models in studies of solute transport in soil is dependent on the accuracy to which the transport parameters involved in this process can be specified. In general, these parameters are determined by using nonlinear optimization techniques to fit theoretical models to experimental data obtained in miscible displacement experiments. The use, therefore, of high-performance software to evaluate these parameters appears to be advantageous since, in addition to consistency and the availability of pre-existing numerical routines, such software allows for the incorporation of new routines that are specifically tailored for the phenomena being simulated. The main objective of this study was to develop a MATLAB 6.5 computational routine to optimize the evaluation of retardation factor (R) and dispersion coefficient (D) from miscible displacement data. This routine was applied to experimental data from three miscible displacement studies of potassium transport in columns filled with an Oxisol of sandy phase, with the setting quality being assessed by the coefficient of accuracy. The proposal routine gave excellent results, reinforcing the consistency of the numerical method used, and indicating that this routine can contribute to the advancement of theoretical studies of the solute-water dynamics in unsaturated porous media.

Keywords: breakthrough curves, computational routine, miscible displacement.

INTRODUÇÃO

O transporte de agroquímicos no solo pela água, para regiões localizadas abaixo do sistema radicular, além de caracterizar prejuízo econômico, pode comprometer a qualidade dos aquíferos subterrâneos. Sendo assim, o conhecimento da dinâmica da água e solutos no solo, em condições de não saturação é importante para a proposição de técnicas e medidas que, além de ganhos econômicos, viabilizem a exploração racional e a preservação dos recursos hídricos (ANAMI et al. 2008; SUÁREZ, 2007).

Estudos vêm buscando descrever, mediante modelos matemáticos, os processos envolvidos no transporte de íons ao longo do perfil do solo (BORGES JUNIOR & FERREIRA, 2006). O sucesso na utilização destes modelos está condicionado à qualidade do ajuste de parâmetros de transporte de solutos. Uma das maneiras de se estimar esses parâmetros mediante experimentos de deslocamento miscível, que resulta em uma distribuição do traçador, que é influenciada pela velocidade microscópica da solução no solo, pelas taxas de difusão e por outras reações químicas (ROSSI et al., 2007; GONÇALVES et al., 2008). De maneira mais geral, o deslocamento miscível pode ser entendido como sendo um processo no qual um fluido se mistura a outro e que se desloca ao longo de um meio poroso. A lixiviação de sais no solo, juntamente com a água que infiltra no solo proveniente da chuva e/ou de irrigação são exemplos de deslocamento miscível (BORGES JUNIOR & FERREIRA, 2006).

Desta forma, diante de estudos de deslocamento miscível, é possível conhecer características de uma determinada substância química, tais como: sua interação com o meio; sua movimentação, e sua persistência no solo. O conhecimento destas características, associado ao uso de modelos de simulação, proporciona melhor entendimento dos processos associados ao deslocamento de solutos no solo, caracterizando-se, portanto, como uma importante ferramenta para a avaliação de riscos de contaminação e dos impactos que podem causar ao meio ambiente.

Desta forma, em trabalhos cujo objetivo é a modelagem matemática e a simulação da dinâmica da água e de solutos no solo, a condução de ensaios experimentais de deslocamento miscível e o posterior ajuste de modelos teóricos às curvas de eluição obtidas é uma prática bastante comum (MIRANDA et al., 2005; MILFONT et al., 2008; RIVERA et al.,2008). Nestes trabalhos, de maneira geral, os dados experimentais são ajustados aos modelos teóricos pela resolução de um problema de otimização não linear em que podem ser determinados, dentre outros, os valores médios da velocidade da água nos poros do solo (v), o coeficiente de dispersão hidrodinâmica (D) e o fator de retardamento (R). Neste sentido, este trabalho teve como objetivo a elaboração de uma rotina na construção de um programa computacional, em ambiente MATLAB 6.5, para otimização não linear no ajuste dos parâmetros de transporte de solutos (v, D e R).

MATERIAL E MÉTODOS

Sob condições de equilíbrio dinâmico, o transporte de solutos no solo é descrito pela eq.(1), cujas unidades são apresentadas em termos dimensionais (M = massa; L = comprimento, e T = tempo):

em que,

C - concentração do soluto na solução do solo, M L^-3;

t - tempo, T;

z - profundidade do solo, L;

R - fator de retardamento (-), e

D - coeficiente de dispersão hidrodinâmica longitudinal, L² T^-1.

A distribuição da concentração de solutos em uma coluna de solo com uma concentração inicial uniforme, C_i, sob equilíbrio dinâmico, a qual será substituída em um tempo "t" igual a zero, por uma solução de concentração constante, C₀, pode ser descrita mediante uma solução analítica da eq.(1) (WRAITH & OR, 1998). Dadas as condições iniciais e de contorno descritas pelas eqs.(2), (3) e (4), a solução analítica da eq.(1) é dada por:

em que,

C(x,t) - concentração do soluto no efluente da coluna em função do tempo t, M L^-3;

C_i - concentração inicial do soluto no solo, M L^-3;

C₀ - concentração do soluto na solução aplicada à coluna, M L^-3;

x - comprimento da coluna, L;

R - fator de retardamento (-), e a denominação "erfc" é a função erro complementar dada por:

Nos casos em que Ci é desprezível, a eq.(5) pode ser simplificada para:

O uso das eqs.(5) e (7) como solução analítica da eq.(1) apresenta resultados muito semelhantes em experimentos de deslocamento miscível (WRAITH & OR, 1998). Assim, na rotina computacional desenvolvida, considerou-se apenas a solução simplificada descrita pela eq.(7):

De acordo com sua própria definição, a função erfc(x) só deve ser calculada nos casos em que . Nos casos em que x < 0, pode-se utilizar a seguinte propriedade:

em que,

O ajuste dos parâmetros de transporte R e D foi feito pelo método dos mínimos quadrados. Neste método numérico, busca-se minimizar a soma dos quadrados dos erros (SQE), pela função objetiva dada por (BURDEN & FAIRES, 2008):

em que,

Cr - concentração relativa medida no efluente da coluna (M L^-3);

- concentração relativa estimada pelo modelo (M L^-3), e

k - número de avaliações da concentração relativa no efluente da coluna ao longo do tempo em que o experimento foi realizado.

Neste sentido, considerando-se a eq.(10), devem-se encontrar os parâmetros R e D, tais que:

Combinando-se as eqs.(7) e (10), tem-se:

A eq.(15) submetida às condições descritas pelas eqs.(11) e (12), gera o seguinte sistema não linear:

A resolução do sistema apresentado pela eq.(16) foi realizada pelo método de Newton. Este método segue o mesmo princípio do método de Newton Raphson, porém com as adaptações necessárias, uma vez que trabalha com conjuntos de funções, de duas ou mais variáveis, para as quais se determina um vetor que seja o "ZERO" das mesmas. Dado um valor inicial (R, D)₀, a cada iteração do método de Newton deve-se resolver o sistema linear descrito pela eq.(17):

As eqs.(18), (19), (20), (21) e (22) descrevem as aproximações obtidas pelo método de diferenças finitas utilizadas no cálculo das derivadas que compõem o sistema linear descrito na eq.(17):

Uma vez estabelecidas as aproximações teóricas, desenvolveu-se, em ambiente MATLAB 6.5, uma rotina computacional que, de acordo com os critérios de parada já discutidos, elabora e resolve o sistema descrito pela eq.(1)7. No ambiente Matlab 6.5, a função erfc foi calculada de acordo com uma rotina interna desenvolvida a partir da metodologia proposta por CODY (1969). Como ocorre em qualquer método iterativo, devem-se estabelecer critérios de parada para se aceitar um ponto (R, D)^k como aproximação para solução (R, D)*. Foram utilizados dois critérios de parada, sendo o primeiro limitando-se em 100 o número de iterações, e o segundo descrito pela eq.(30) dada por:

A convergência teórica do método de Newton é uma convergência local, ou seja, para que exista, a aproximação inicial escolhida deve estar suficientemente próxima da solução exata do sistema não linear que se queira estudar (BURDEN & FAIRES, 2008). Neste sentido, criou-se como primeiro passo, no programa desenvolvido, um algoritmo que, a partir de diferentes combinações de R e D, aceita como "aproximação inicial" aquele par associado ao menor SQE (Figura 1).

A aproximação inicial, (R, D)⁰, estabelecida pelo fluxograma apresentado na Figura 1, foi utilizada como dado de entrada no programa principal, no qual, com base no fluxograma apresentado na Figura 2, realiza-se o ajuste por mínimos quadrados dos parâmetros de transporte R e D.

A rotina computacional proposta foi utilizada para ajustar os parâmetros de transporte aos dados experimentais de SOUZA (2007) e de GONÇALVES (2007). Ambos os autores realizaram três repetições de ensaios de curvas de eluição do íon potássio em colunas preenchidas com um Latossolo Vermelho-Amarelo, textura média, classificado de acordo com o Sistema Brasileiro de Classificação de solos (EMBRAPA, 1999).

Em seu estudo, SOUZA (2007) conduziu a aplicação de uma solução de 1.500 ppm de potássio, mediante uma fonte de nitrato de potássio (KNO₃). GONÇALVES (2007), por sua vez conduziu a aplicação de uma solução de 2.000 ppm de potássio, mediante a combinação de dois sais: o nitrato de potássio (KNO₃) e o cloreto de potássio (KCl). Os ajustes obtidos pela rotina computacional desenvolvida no Matlab 6.5 foram comparados com os ajustes obtidos pelo programa Disp (BORGES JUNIOR & FERREIRA, 2006). A qualidade dos ajustes obtidos foi avaliada mediante a utilização do índice de desempenho (c) proposto por CAMARGO & SENTELHAS (1997) (eq.(24).

em que,

r - coeficientes de correlação de Pearson, e

d - coeficiente de concordância.

O coeficiente de concordância (WILLMOTT et al., 1985) é descrito por:

em que,

n - número de pontos avaliados experimentalmente no ensaio de deslocamento miscível;

Cr - concentração relativa obtida experimentalmente;

Cr^*- concentração relativa estimada, e

- média das concentrações relativas obtidas experimentalmente.

O desempenho dos métodos, avaliados por meio do índice (c), foi classificado como ótimo (O), muito bom (MB), bom (B), mediano (Md), sofrível (S), mau (M) e péssimo (P), Tabela1.

RESULTADOS E DISCUSSÃO

Para cada um dos seis ensaios experimentais considerados, foram obtidos os valores dos parâmetros de transporte ajustados, utilizando-se do Matlab 6.5 e do programa Disp (BORGES JUNIOR & FERREIRA, 2006) (Tabela 2).

Pôde-se verificar de maneira geral que, para um mesmo ensaio experimental, existem diferenças entre os parâmetros ajustados pelas duas rotinas computacionais utilizadas. Tais diferenças justificam-se pelas metodologias empregadas no calculo da função erro complementar (erfc) e pela simplificação da eq.(5), bem como, pelo fato de que SOUZA (2007) aplicou uma solução de 1.500 ppm de potássio, via fonte de nitrato de potássio (KNO³), e GONÇALVES (2007) aplicou uma solução de 2.000 ppm de potássio, mediante a combinação de dois sais: o nitrato de potássio (KNO³) e o cloreto de potássio (KCl). Dessa forma, o efeito do íon acompanhante, no caso de GONÇALVES (2007), por se tratar de dois elementos com carga negativa (Nitrato e Cloreto), foi evidenciado maior efeito de adsorção do íon potássio no solo. Portanto, os parâmetros evidenciaram, numericamente, o efeito físico do processo.

Assim sendo, mesmo quando as informações apresentadas pela Tabela 2 são utilizadas no ajuste das curvas de distribuição de efluentes aos dados experimentais de SOUZA (2007) e aos dados experimentais de GONÇALVES (2007), constatou-se que a diferença numérica dos parâmetros ajustados tem pouca influência no comportamento do processo modelado, mas que evidencia fisicamente o processo da dinâmica do íon potássio no solo.

Mediante as Figuras 3A, 3B e 3C (referentes ao ensaio de SOUZA, 2007) e às Figuras 3D, 3E e 3F (referentes ao ensaio de GONÇALVES, 2007), no que se refere à qualidade dos ajustes, pôde-se constatar, visualmente, que os ajustes do modelo teórico mostraram boa aderência aos dados experimentas. Assim também, se forem analisados os valores do parâmetro fator de retardamento (considerando o valor de C/Co a 0,5), pode-se verificar que tais valores são semelhantes aos valores dos ajustes obtidos para o parâmetro fator de retardamento (R) (Tabela 2), evidenciando a qualidade do ajuste, tanto numérica, quanto na representação física do processo de dinâmica do íon potássio.

 

Tal observação pode ser confirmada pela análise da Tabela 3, na qual cada rotina computacional foi avaliada, nos diferentes ensaios experimentais considerados, e foram apresentados os coeficientes de correlação (r), de concordância (d) e o de exatidão (c). Verificou-se que, em todos os casos avaliados, ambos os modelos tiveram seu desempenho classificado como "ótimo", conforme critério de CAMARGO & SENTELHAS (1997).

CONCLUSÕES

O modelo proposto neste trabalho para otimização de parâmetros apresentou um desempenho ótimo, o que reforça a consistência do método numérico utilizado.

Recebido pelo Conselho Editorial em: 10-3-2010

Aprovado pelo Conselho Editorial em: 11-7-2011

Referências bibliográficas

• ANAMI, M.H.; SAMPAIO, S.C.; SUSZEK, M.; GOMES, S.D.; QUEIROZ, M.M.F. Deslocamento miscível de nitrato e fosfato proveniente de água residuária da suinocultura em colunas de solo. Revista Brasileira de Engenharia Agrícola e Ambiental, Campina Grande, v.12, n.1, fev. 2008.
• BORGES JUNIOR, J.C.F.; FERREIRA, P.A. Equações e programa computacional para cálculo do transporte de solutos no solo. Revista Brasileira de Engenharia Agrícola e Ambiental, Campina Grande, v.10, n.3, p.604-611, 2006.
• BURDEN, R.L.; FAIRES, J.D. Análise numérica 8.ed. São Paulo: Cengage Learning, 2008. 736 p.
• CAMARGO, A.P.; SENTELHAS, P.C. avaliação do desempenho de diferentes métodos de estimativa da evapotranspiração potencial no estado de São Paulo, Brasil. Revista Brasileira de Agrometeorologia, Santa Maria, v.5, n.1, p.89-97,1997.
• CODY, W.J. Rational Chebyshev Approximations for the Error Function. Mathematics Computation, Providence, p.631-638, 1969.
• EMBRAPA. EMPRESA BRASILEIRA DE PESQUISA AGROPECUÁRIA. Sistema brasileiro de classificação dos solos Rio de janeiro: Embrapa-SPI, 1999. 412 p.
• GONÇALVES, A.D.M.A. Efeito da temperatura no transporte dos íons potássio e nitrato no solo 2007. 82 f. (Dissertação de Mestrado) - Escola Superior de Agricultura "Luiz de Queiroz", Universidade de São Paulo, Piracicaba, 2007.
• GONÇALVES, A.D.M.A.; MIRANDA, J.H.; ROSSI, P.; SABADIN, J.F.G.; KAMOGAWA, M.Y. Temperature effect in potassium and nitrate ions in soil transport. Engenharia Agrícola, Jaboticabal, v.28, n.3, p.438-447, set. 2008.
• MIRANDA, J.H.; DUARTE, S.N.; LIBARDI, P.L.; FOLEGATTI, M.V. Simulação do deslocamento de potássio em colunas verticais de solo não saturado. Engenharia Agrícola, Jaboticabal, v.25, n.3, p.677-685, dez. 2005.
• MILFONT, M.L.; ANTONINO, A.C.D.; MARTINS, J.M.F.; NETTO, A.M.; GOUVEIA, E.R.; CORREA, M.M. Transporte do paclobutrazol em colunas de solos. Revista Brasileira de Ciência do Solo, Viçosa-MG, v.32, n.5, out. 2008.
• RIVERA, R.N.C.; MIRANDA, J.H.; DUARTE, S.N.; BOTREL, T.A. Modelo aplicado à dinâmica da água e do potássio no solo sob irrigação por gotejamento: análise de sensibilidade. Engenharia Agrícola, Jaboticabal, v.28, n.3, set. 2008.
• ROSSI, P.; MIRANDA, J.H.; DUARTE, S.N. Curvas de distribuição de efluentes do íon nitrato em amostras de solo deformadas e indeformadas. Engenharia Agrícola, Jaboticabal, v.27, n.3, dez. 2007.
• SOUZA, A.L. Teoria de conjuntos fuzzy no estudo da dinâmica da água e de solutos no solo 2007. 157 f. Tese (Doutorado) - Faculdade de Engenharia Agrícola, Universidade Estadual de Campinas, Campinas, 2007.
• SUÁREZ, F.; BACHMANN, J.; MUÑOZ, J.F.; ORTIZ, C.; TYLER, S.W.; ALISTER, C.; KOGAN, M. Transport of simazine in unsaturated sandy soil and predictions of its leaching under hypothetical field conditions. Journal of Contaminant Hydrology, Exeter, v.94, n.3-4, p.166-177, 2007.
• WILLMOTT, C.J.; ACKLESON, S.G.; DAVIS, R.E.; FEDDEMA, J.J.; KLINK, K.M.; LEGATES, D.R.; O'DONNELL, J.; ROWE, C.M. Statistics for the evaluation and comparision of models. Journal of Geophysical Research, Ottawa, v. 90, n.C5, p.8.995-9.005, 1985.
• WRAITH, J.M.; OR, D. Nonlinear parameter estimation using spreadsheet software. Journal of Natural Resources, Albuquerque, v.27, p.13-19, 1998.

Datas de Publicação

• Publicação nesta coleção
  09 Jan 2012
• Data do Fascículo
  Dez 2011

Histórico

• Aceito
  11 Jul 2011
• Recebido
  10 Mar 2010