2018-08-19 18:47:08 +00:00
|
|
|
#!/usr/bin/env python3
|
|
|
|
import sys
|
|
|
|
import os
|
|
|
|
import numpy as np
|
|
|
|
Nx = int(sys.argv[1])
|
|
|
|
Ny = int(sys.argv[2])
|
|
|
|
param = sys.argv[3]
|
|
|
|
|
|
|
|
buffer = bytearray(b" "*(3*Nx*Ny))
|
|
|
|
|
|
|
|
|
|
|
|
Lx = 32 #128 #Physikalische Laenge des Gebietes in x-Richtung
|
|
|
|
Ly =32 #128#Physikalische Laenge des Gebietes in y-Richtung
|
|
|
|
|
|
|
|
x, y = np.meshgrid(np.arange(Nx)* Lx/Nx, np.arange(Ny)* Ly/Ny) #x-Array
|
|
|
|
kx, ky = np.meshgrid(np.fft.fftfreq(Nx,Lx/(Nx*2.0*np.pi)), np.fft.fftfreq(Ny,Ly/(Ny*2.0*np.pi)))
|
|
|
|
ksq = kx*kx + ky*ky
|
|
|
|
|
|
|
|
c = 0.0+(np.random.random((Ny,Nx))-0.5)*0.1
|
2018-08-22 17:43:16 +00:00
|
|
|
eps=0.3
|
|
|
|
delta=0.0
|
|
|
|
#eps = 0.1
|
|
|
|
#delta = 1.0
|
2018-08-19 18:47:08 +00:00
|
|
|
|
|
|
|
t = 0.0
|
|
|
|
dt = 0.01
|
|
|
|
T_End = 300000
|
|
|
|
N_t = int(T_End / dt)
|
|
|
|
|
|
|
|
plotEveryNth = 100
|
|
|
|
ck = np.fft.fft2(c) #FFT(c)
|
|
|
|
# Lineare Terme
|
|
|
|
def rhs_lin(ksq):
|
|
|
|
result=eps-(1.0-ksq)**2
|
|
|
|
return result
|
|
|
|
|
|
|
|
Eu=1.0/(1.0-dt*rhs_lin(ksq))
|
|
|
|
i = 0
|
|
|
|
def lerp_sat(a, t, b):
|
|
|
|
v = (1-t)*a+t*b
|
|
|
|
v = int(v)
|
|
|
|
if v < 0:
|
|
|
|
v = 0
|
|
|
|
if v > 255:
|
|
|
|
v = 255
|
|
|
|
return v
|
|
|
|
|
|
|
|
while True:
|
|
|
|
i+= 1
|
|
|
|
ck=Eu*(ck+dt*(delta*np.fft.fft2(np.fft.ifft2(ck)**2)-np.fft.fft2(np.fft.ifft2(ck)**3)))
|
|
|
|
c=np.fft.ifft2(ck)
|
2018-08-22 17:43:16 +00:00
|
|
|
eps = 0.1+0.2*np.cos(i/10000)
|
|
|
|
delta = np.sin(i/10000)
|
2018-08-19 18:47:08 +00:00
|
|
|
if(i % plotEveryNth == 0):
|
|
|
|
myc = c.real
|
|
|
|
myc = (myc-np.min(myc))/(np.max(myc)-np.min(myc))
|
|
|
|
for px in range(Nx):
|
|
|
|
for py in range(Ny):
|
|
|
|
idx = 3*(px+Nx*py)
|
|
|
|
buffer[idx+0] = lerp_sat(0xff, myc[py,px], 0x00)
|
|
|
|
buffer[idx+1] = lerp_sat(0x00, myc[py,px], 0xff)
|
|
|
|
buffer[idx+2] = 0
|
|
|
|
os.write(1, buffer)
|