Exibir código-fonte para Equação de Cahn-Hilliard

'''Grupo: Arthur Dornelles, Bruno Zanette, Gabriel De David, Guilherme Hoss'''

O objetivo deste trabalho é resolver computacionalmente a equação de Cahn-Hilliard, que descreve o processo de decomposição spinodal de uma mistura binária, utilizando o método FTCS (''Forward Time Centered Space'').

== Decomposição Espinodal ==

Decomposição espinodal é o nome dado ao processo no qual uma pequena perturbação de um sistema faz com que, uma fase homogênea termodinamicamente instável, diminua sua energia e separe-se espontaneamente em duas outras fases coexistentes, esse é um processo que ocorre sem nucleação, ou seja, é instantâneo. Ela é observada, por exemplo, em misturas de metais ou polímeros e pode ser modelada pela equação de Cahn-Hilliard.

==A Equação de Cahn-Hilliard==

A equação de Cahn-Hilliard descreve o processo de decomposição espinodal de uma mistura binária.
Em outras palavras, é uma equação que descreve o processo de separação de fase entre dois componentes de um fluido binário que se separam de maneira espontânea.

Consideraremos - de início - uma mistura binária de dois componentes A e B descritas pelas concentrações dos fluidos <math> c_a(x,t) </math> e <math> c_b(x,t) </math>, respectivamente.

Além disso, podemos considerar que - para uma mistura binária - <math> c_a(x,t) + c_b(x,t) = 1</math> e portanto podemos simplificar para apenas uma concentração <math> c(x,t) </math>:
:<math>
c_a(x,t) = c(x,t), c_b(x,t) = 1 - c(x,t)
</math>

Tendo isso em vista, podemos agora utilizar a primeira lei de Fick da difusão:
:<math>
J = -D\nabla c
</math>

juntamente da equação da continuidade:
:<math>
\frac{\partial c}{\partial t} + \nabla \cdot \vec J = 0
</math>

Onde <math>D</math> é o coeficiente de difusão e <math>\vec J</math> é o fluxo de difusão. Em seguida, ao combinarmos ambas as equações anteriores o resultado gera a segunda lei de Fick da difusão:

:<math>
\frac{\partial c}{\partial t} = D {\nabla}^2 c
</math>

Por definição, verificou-se que a concentração não poderia ser a razão da difusão, portanto outra força estaria presente. E, nesse caso, encontrou-se que a principal força responsável pela difusão negativa é o potencial químico. Portanto, outra equação pode ser derivada para generalizar a primeira lei de Fick:

:<math>
J = -M \nabla \mu
</math>

Onde <math>M</math> é a mobilidade das partículas (análoga à D) e <math>\mu</math> é o potencial químico.
Com essa nova equação podemos agora também deduzir uma nova equação para a segunda lei de Fick:

:<math>
\frac{\partial c}{\partial t} = M {\nabla}^2 \mu
</math>

Essa equação também é conhecida como equação de Cahn-Hilliard

Nessa equação, podemos usar a definição do potencial químico através da densidade de energia livre de Gibbs como:

:<math>
\mu = \frac{\partial g}{\partial c}
</math>

Onde <math>g</math> é a densidade de energia livre de Gibbs e <math>c</math> é a concentração.

Tendo em vista a substituição do termo <math>\mu</math> por um termo que envolva a concentração dos fluidos, utiliza-se uma equação que descreve a densidade de energia desse sistema através da concentração dos mesmos:

:<math>
F = \int_{V}^{} f(c) + {\kappa |\nabla c|}^2 dV
</math>

Nesse caso, <math>F</math> é a energia livre de Gibbs, <math>f(c)</math> é a densidade de energia livre devido à contribuições de ambas as fases homogêneas e <math>{\kappa|\nabla c|}^2</math> é a densidade de energia livre devido ao gradiente de concentração na interface (ou energia de interface).

Além disso, a função <math>f(c)</math> tem o formato de um poço de potencial duplo, que pode ser representado pela seguinte equação:
:<math>
f(c) = \frac{(c^2 - 1)^2}{4}
</math>

Levando essas informações em conta e - utilizando um parâmetro <math>\epsilon</math> análogo à largura da interface - que é descrito por <math>\kappa = \epsilon ^2</math> é possível encontrar uma equação que descreve a densidade de energia lvire de Gibbs para um sistema duplo-fásico:
:<math>
g(c) = f(c) + {\epsilon}^2{|\nabla c |}^2
</math>

Com essas igualdades agora se torna possível o cálculo de <math>\mu</math> em função da concentração dos fluidos:

:<math>
\mu = \frac{\partial g}{\partial c} = \frac{\partial f(c)}{\partial c} + \frac{\partial {(\epsilon}^2{|\nabla c |}^2)}{\partial c} = \frac{\partial \frac{((c^2 - 1)^2}{4})}{\partial c} + {\epsilon}^2 {\nabla}^2 c = c^3 - c + {\epsilon}^2 {\nabla}^2 c
</math>

== Método FTCS (Forward Time Centered Space) ==

O FTCS é um método numérico utilizado para resolver equações diferenciais parciais, tais como a difusão do calor e do transporte de massa, traduzindo, significa "Progressivo no tempo, avançado no espaço". Esse método pode ser utilizado em sua forma implícita ou explícita que estão descritas abaixo.

:<math>n\to\Delta t</math>

:<math>j\to\Delta x</math>

===FTCS Explicito===

:<math>\frac{\partial f}{\partial t}\to \frac{f_{j}^{n+1}-f_{j}^{n}}{\Delta t}</math>

:<math>\frac{\partial ^2 f}{\partial x^2}\to \frac{f_{j-1}^{n}-2 f_{j}^{n} + f_{j+1}^n}{\Delta x^2}</math>

Para difusão:

:<math>f_j^{n+1}= f_j^{n} + \frac{D\Delta t}{(\Delta x)^2} (f_{j-1}^{n} - 2f_j^{n} + f_{j+1}^{n})</math>

=== FTCS Implicito (BTCS) ===

:<math>\frac{\partial f}{\partial t}\to \frac{f_{j}^{n}-f_{j}^{n-1}}{\Delta t}</math>

:<math>\frac{\partial ^2 f}{\partial x^2}\to \frac{f_{j-1}^{n}-2 f_{j}^{n} + f_{j+1}^n}{\Delta x^2}</math>


Para difusão:

:<math>f_j^{n+1}= f_j^{n} + \frac{D\Delta t}{(\Delta x)^2} (f_{j-1}^{n+1} - 2f_j^{n+1} + f_{j+1}^{n+1})</math>


== Resolução do Cahn-Hilliard Equation para FTCS explicito para x somente: ==

:<math>\displaystyle \frac{\partial c}{\partial t} = D\nabla^{2}(c^{3}-c-\gamma^2\nabla^{2}c)</math>  

:<math>\displaystyle \frac{c_{j}^{n+1}-c_{j}^{n}}{\Delta t} =  D\displaystyle \frac{\partial ^2 }{\partial x^2}(c^3 - c - \gamma^2 \displaystyle \frac{\partial ^2 c}{\partial x^2})</math>

:<math>\displaystyle \frac{c_{j}^{n+1}-c_{j}^{n}}{\Delta t} = D\frac{u_{j-1}^n-2u_j^n + u_{j+1}^n}{(\Delta x)^2}</math>


:<math>\displaystyle \frac{c_{j}^{n+1}-c_{j}^{n}}{\Delta t} = D\left(\frac{(c_{j-1}^n)^3-2(c_j^n)^3 + (c_{j+1}^n)^3}{(\Delta x)^2} - \frac{c_{j-1}^n-2 c_i^n + c_{j+1}^n}{(\Delta x)^2} - \gamma^2\frac{c_{j-2}^n-4c_{j-1}^n + 6c_{j}^n -4c_{j+1}^n + c_{j+2}^n}{(\Delta x)^4}\right)</math>

:<math> c_{j}^{n+1} = D\Delta t \left (\frac{(c_{j-1}^n)^3-2(c_i^n)^3 + (c_{j+1}^n)^3}{(\Delta x)^2} - \frac{c_{j-1}^n-2 c_j^n + c_{j+1}^n}{(\Delta x)^2} - \gamma^2\frac{c_{j-2}^n-4c_{j-1}^n + 6c_{j}^n -4c_{j+1}^n + c_{j+2}^n}{(\Delta x)^4} \right) + c_j^n</math>

:<math>c_{j}^{n+1} = \frac{D\Delta t}{(\Delta x)^2} \left ((c_{j-1}^n)^3-2(c_i^n)^3 + (c_{j+1}^n)^3 - {c_{j-1}^n+2 c_j^n - c_{j+1}^n} - \gamma^2\frac{c_{j-2}^n-4c_{j-1}^n + 6c_{j}^n -4c_{j+1}^n + c_{j+2}^n}{(\Delta x)^2} \right) + c_j^n </math>

'''Condição de Estabilidade '''

:<math>\Delta t < \displaystyle\frac{(\Delta x)^2}{2}</math>

== Referências ==