Грей Скотт Модель - PullRequest
Новая
       9

Грей Скотт Модель

1 голос
/ 21 апреля 2019

Я пытаюсь показать модель диффузии Грея Скотта. Я продолжаю получать кучу предупреждений об ошибках во время выполнения, хотя чувствую, что мой код действительно близок к исправлению. Что-то не так с моей дискретизацией? 

import numpy as np
import matplotlib.pyplot as plt

#parameters
N=128
F=.042
k=.062
Du=(2**-5)*(N**2/6.25)
Dv=(1**-5)*(N**2/6.25)
tend=1000                             
dt=tend/N
t=0

#start arrays
U=np.ones((N,N))
V=np.zeros((N,N))

#Initial Value Boxes (20x20 in middle)
low=int(((N/2)-10))
high=int(((N/2)+10))+1
U[low:high,low:high]=.5
V[low:high,low:high]=.25

#Random Noise
U+=.01*np.random.random((N,N))
V+=.01*np.random.random((N,N))

#Laplace
def Laplace(f):
    return np.roll(f,1)+np.roll(f,-1)+np.roll(f,1,axis=False)+np.roll(f,-1,axis=False)-4*f

#Solve
pstep=100
for t in range(tend):
    U+=((Du*Laplace(U))-(U*V*V)+(F*(1-U)))
    V+=((Dv*Laplace(V))+(U*V*V)-(V*(F+k)))
    if t%pstep ==0:
        print(t//pstep)
        plt.imshow(U, interpolation='bicubic',cmap=plt.cm.jet)

1 Ответ

1 голос
/ 22 апреля 2019

Хорошо, я заставил его работать, изменив несколько вещей в расчете, но в основном изменив числовую стабильность, значительно уменьшив коэффициенты диффузии и уменьшив временной шаг. Конечным результатом этого является то, что все моделирование меняется меньше между каждым шагом, поэтому значение изменения намного меньше.

Ошибки, которые вы получали, были связаны с переполнением чисел в вычислениях dU и dV, поскольку при замедлении всего (больше временных шагов) вам не нужны такие массивные числа в dU и dV

import numpy as np
import matplotlib.pyplot as plt

# parameters
N = 128
F = .042
k = .062
# Du=(2**-5)*(N**2/6.25) # These were way too high for the 
# numeric stability given the timestep
Du = 0.1
# Dv=(1**-5)*(N**2/6.25)
Dv = 0.5
tend = 1000
dt = tend / N
t = 0
dt = 0.1  # Timestep - 
# the smaller you go here, the faster you can let the diffusion coefficients be

# start arrays
U = np.ones((N, N))
V = np.zeros((N, N))

# Initial Value Boxes (20x20 in middle)
low = int(((N / 2) - 10))
high = int(((N / 2) + 10)) + 1
U[low:high, low:high] = .5
V[low:high, low:high] = .25

# Random Noise
U += .01 * np.random.random((N, N))
V += .01 * np.random.random((N, N))


# Laplace
def Laplace(f):
    return np.roll(f, 1) + np.roll(f, -1) + np.roll(f, 1, axis=False) + np.roll(f,-1,                                                                              axis=False) - 4 * f


# Solve
pstep = 100
for t in range(tend):
    dU = ((Du * Laplace(U)) - (U * V * V) + (F * (1 - U))) * dt
    dV = ((Dv * Laplace(V)) + (U * V * V) - (V * (F + k))) * dt
    U += dU
    V += dV
    if t % pstep == 0:
        print(t // pstep)
        plt.imshow(U, interpolation='bicubic', cmap=plt.cm.jet, vmin=0, vmax=1)

Конечно, изменения, которые я сделал, немного изменяют физику, и вам нужно будет изменить свои t и pstep, чтобы эти значения имели смысл. А также проверьте, как вы рассчитывали Du и Dv. Если эти значения на самом деле должны быть ~ 8000, то вам нужен намного меньший временной шаг.

Для справки, объяснение модели Грея Скотта здесь

...