1# ##### BEGIN GPL LICENSE BLOCK ##### 2# 3# This program is free software; you can redistribute it and/or 4# modify it under the terms of the GNU General Public License 5# as published by the Free Software Foundation; either version 2 6# of the License, or (at your option) any later version. 7# 8# This program is distributed in the hope that it will be useful, 9# but WITHOUT ANY WARRANTY; without even the implied warranty of 10# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 11# GNU General Public License for more details. 12# 13# You should have received a copy of the GNU General Public License 14# along with this program; if not, write to the Free Software Foundation, 15# Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. 16# 17# ##### END GPL LICENSE BLOCK ##### 18 19import numpy as np 20try: 21 from numba import jit 22 23 @jit 24 def numba_reaction_diffusion(n_verts, n_edges, edge_verts, a, b, diff_a, diff_b, f, k, dt, time_steps): 25 arr = np.arange(n_edges)*2 26 id0 = edge_verts[arr] # first vertex indices for each edge 27 id1 = edge_verts[arr+1] # second vertex indices for each edge 28 for i in range(time_steps): 29 lap_a = np.zeros(n_verts) 30 lap_b = np.zeros(n_verts) 31 lap_a0 = a[id1] - a[id0] # laplacian increment for first vertex of each edge 32 lap_b0 = b[id1] - b[id0] # laplacian increment for first vertex of each edge 33 34 for i, j, la0, lb0 in zip(id0,id1,lap_a0,lap_b0): 35 lap_a[i] += la0 36 lap_b[i] += lb0 37 lap_a[j] -= la0 38 lap_b[j] -= lb0 39 ab2 = a*b**2 40 #a += eval("(diff_a*lap_a - ab2 + f*(1-a))*dt") 41 #b += eval("(diff_b*lap_b + ab2 - (k+f)*b)*dt") 42 a += (diff_a*lap_a - ab2 + f*(1-a))*dt 43 b += (diff_b*lap_b + ab2 - (k+f)*b)*dt 44 return a, b 45except: 46 pass 47