Simulaçao Belousov-Zhabotinsky: mudanças entre as edições

De Física Computacional
Ir para navegação Ir para pesquisar
Sem resumo de edição
Sem resumo de edição
 
(Uma revisão intermediária por um outro usuário não está sendo mostrada)
Linha 22: Linha 22:
def inicia():
def inicia():
     global u, v, nextu, nextv
     global u, v, nextu, nextv
     u = np.loadtxt('u_file.txt')
     u = np.loadtxt('u_file.txt') ## condições iniciais
     v = np.loadtxt('v_file.txt')
     v = np.loadtxt('v_file.txt') ## condições iniciais


     nextu = zeros([n, n])
     nextu = zeros([n, n])
Linha 63: Linha 63:
     ax.axis(False);
     ax.axis(False);


# In[ ]:





Edição atual tal como às 14h32min de 30 de março de 2021

<source lang = "py"> import matplotlib matplotlib.use('TkAgg') from pylab import * import numpy as np

n = 100 Dh = 1./n Dt = 0.001

Du = 1e-5 Dv = 1e-5

eps = 0.2

q = 1e-3

f = 1


def inicia():

   global u, v, nextu, nextv
   u = np.loadtxt('u_file.txt') ## condições iniciais
   v = np.loadtxt('v_file.txt') ## condições iniciais
   nextu = zeros([n, n])
   nextv = zeros([n, n])

def update():

   global u, v, nextu, nextv
   for i in range(1, n - 1):
       for j in range(1, n - 1):
           uC = u[i, j]
           uR = u[(i + 1) % n, j]
           uL = u[(i - 1) % n, j]
           uU = u[i, (j + 1) % n]
           uD = u[i, (j - 1) % n]
           vC = v[i, j]
           vR = v[(i + 1) % n, j]
           vL = v[(i - 1) % n, j]
           vU = v[i, (j + 1) % n]
           vD = v[i, (j - 1) % n]
           uLap = (uR + uL + uU + uD - 4 * uC) / (Dh ** 2)
           vLap = (vR + vL + vU + vD - 4 * vC) / (Dh ** 2)
           nextu[i, j] = uC + (uC * (1 - uC) - f * vC * (uC - q) / (uC + q) + Du * uLap) * Dt / eps
           nextv[i, j] = vC + (uC - vC + Dv * vLap) * Dt
   u, nextu = nextu, u
   v, nextv = nextv, v

def observeu():

   global u, v, nextu, nextv
   fig, ax = plt.subplots(figsize=(16, 9))
   ax.matshow(u, cmap='plasma')
   ax.axis(False);

def observev():

   global u, v, nextu, nextv
   fig, ax = plt.subplots(figsize=(16, 9))
   ax.matshow(v, cmap='winter')
   ax.axis(False);


inicia() t = 0 observeu() for t in range(1, 100001):

   update()
   if t % 500 == 0:
       print(t)
   if t % 1000 == 0:
       observeu()
       savefig('u-BZ-{}.png'.format(t))
       observev()
       savefig('v-BZ-{}.png'.format(t))