Skip to content

Advection-Diffusion Eq, Turing instabilities linear and non-linear case

Notifications You must be signed in to change notification settings

96mat/Turing-Instability

Folders and files

NameName
Last commit message
Last commit date

Latest commit

 

History

54 Commits
 
 
 
 
 
 

Repository files navigation

Turing Instability

la repo è costituita da due script principalmente:

  • Il codice MixedWF_CH.ipynb, MixedWF_CH_v1.ipynb e MixedWF_CH_v2.ipynb in cui viene approcciato il problema non lineare $(1)$, con risultati non adatti
  • Il codice TN.ipynb in cui viene approcciato un problema simile ma lineare $(2)$, i risultati sembrano ok
  • Plot delle concentrazioni $u,v$ in FreeFEM

MixedWF_CH.ipynb, MixedWF_CH_v1.ipynb, MixedWF_CH_v2.ipynb

Advection-Diffusion Eq

$$\begin{cases} u_t-\varDelta u=f(u,v)\\ v_t -d\varDelta v =g(u,v) \qquad(1)\\ u(x,y,z)=u_0(x,y)\\ v(x,y,z)=v_0(x,y)\\ + \text{Neuman BC}\\ +\text{Periodic BC} \end{cases}$$

dove i termini non lineari sono:

$$\begin{align} f(u,v)&=u-\alpha v +\gamma u v - u^3\\ g(u,v)&= u-\beta v \end{align}$$

and $u_0(x,y),v_0(x,y)$ sono le initial condition e sono inizializzate con la classe IC cioè distribuzioni pseudo-randomiche.

#Initial condtions (Random)
class IC(UserExpression):
    def eval(self,values,x):
        values[0] = 0.1*random() -0.1*random()
        values[1] = 0.1*random() -0.1*random()
    def value_shape(self):
        return(2,)

$\beta, \thinspace d, \thinspace \alpha,\thinspace \gamma$ sono costanti, il risultato:

Ulteriormente

Idea per l'implementazione delle funzioni ff, gg in modo simbolico da poter usare nella formulazione debole:

#Define the forcing terms
u_n1 = variable(u_n1)
u_n2 = variable(u_n2)
ff   = u_n1-alpha*u_n2+gamma*u_n1*u_n2-u_n1**3
gg   = u_n1-beta*u_n2

in alternativa all'interno della formulazione debole, si accedono alle componenti di

u_n=Function(V)
u_n1, u_n2 =split(u_n)

usando u_n[0], u_n[1] e scrivendo direttamente $f(u,v),\thinspace g(u,v)$ nella formulazione debole:

F =  (u_1*v_1/dt)*dx \
    + (u_2*v_2/dt)*dx \
    + inner(grad(u_1), grad(v_1))*dx \
    + d*inner(grad(u_2), grad(v_2))*dx \
    -(dot(u_n1,v_1)/dt + dot(u_n2,v_2)/dt)*dx \
    -(dot(u_n[0]-alpha*u_n[1]+gamma*u_n[0]*u_n[1]-u_n[0]**3,v_1) + dot(u_n[0]-beta*u_n[1],v_2))*dx 

altrimenti

def ff(u_n1,u_n2):
    return u_n1-A*u_n2+G*u_n1*u_n2-u_n1**3
def gg(u_n1,u_n2):
    return u_n1-B*u_n2

# Define variational problem
F =  (u_1*v_1/dt)*dx \
    + (u_2*v_2/dt)*dx \
    + inner(grad(u_1), grad(v_1))*dx \
    + d*inner(grad(u_2), grad(v_2))*dx \
    -(dot(u_n1,v_1)/dt + dot(u_n2,v_2)/dt)*dx \
    -(dot(ff(u_n1,u_n2),v_1) + dot(gg(u_n1,u_n2),v_2))*dx 

MixedWF_CH_v2.ipynb

E' risolto in 2 possibili modi, usando un solver non-lineare alternativo a Newton, ovvero SNES i quali i sistemi lineari sono risolti con GMRES+HYPRE_AMG (preconditioner)

TN.ipynb

$$\begin{cases} \frac{\partial X}{\partial t}=a(X-h)+b(Y-k)+\mu \nabla^2 X \qquad (2)\\\ \frac{\partial Y}{\partial t}=c(X-h)+d(Y-k)+\nu \nabla^2 Y\\ +\text{Periodic BC} \end{cases}$$

la cui Forma debole (pagina 9)

$$\begin{align} \int_{\Omega} (X^{n+1}-X^n)v+\mu \Delta t \nabla X^{n+1}\nabla v- \Delta t(aX^{n+1}+bY^{n+1})v\thinspace d\Omega &=0\\\ \int_{\Omega} (Y^{n+1}-Y^n)v+\nu \Delta t \nabla Y^{n+1}\nabla v- \Delta t(cX^{n+1}+dY^{n+1})v\thinspace d\Omega &=0 \end{align}$$

la soluzione che ho ottenuto:

Ulteriormente

FreeFEM

plot ottenuto attraverso il codice FreeFEM

Ulteriormente

About

Advection-Diffusion Eq, Turing instabilities linear and non-linear case

Resources

Stars

Watchers

Forks

Releases

No releases published

Packages

No packages published