62 lines
1.5 KiB
Python
62 lines
1.5 KiB
Python
|
#!/usr/bin/env python3
|
||
|
import numpy as np
|
||
|
import time
|
||
|
import os
|
||
|
import sys
|
||
|
import colorsys
|
||
|
Nx = int(sys.argv[1])
|
||
|
Ny = int(sys.argv[2])
|
||
|
N = max(Nx, Ny)
|
||
|
|
||
|
def CpxRandNormal(min, max, N):
|
||
|
mu = (min+max)/2
|
||
|
sigma = (max-min)/2
|
||
|
return np.random.normal(mu, sigma, (N, N))+1.0j*np.random.normal(mu, sigma, (N, N))
|
||
|
def CpxRand(min, max, N):
|
||
|
return np.random.uniform(min, max, (N, N))+1.0j*np.random.uniform(min, max, (N, N))
|
||
|
|
||
|
A = CpxRandNormal(-100, 100, N)
|
||
|
B = CpxRand(-1, 1, N)*0.1
|
||
|
C = CpxRand(-1, 1, N)
|
||
|
|
||
|
A_ex = np.linalg.inv(np.eye(N)-B)@C
|
||
|
|
||
|
dir = -1
|
||
|
solve = True
|
||
|
t = 0
|
||
|
while True:
|
||
|
t += 0.003
|
||
|
if solve:
|
||
|
A = 0.9*A+0.1*(B@A+C)
|
||
|
else:
|
||
|
A += CpxRandNormal(0, 0.008, N)
|
||
|
diff = A-A_ex
|
||
|
#angles = (np.angle(diff)+np.pi)/(2*np.pi
|
||
|
#diff = np.abs(diff)
|
||
|
#diff = 1-diff/(1+diff)
|
||
|
a = np.abs(np.real(diff))
|
||
|
a = a/(a+1)
|
||
|
b = np.abs(np.imag(diff))
|
||
|
b = b/(b+1)
|
||
|
img = np.zeros((Ny, Nx, 3))
|
||
|
for x in range(Nx):
|
||
|
for y in range(Ny):
|
||
|
c = colorsys.hsv_to_rgb((t+0.2+a[x, y]*0.2)%1, 1, 0.5*b[x, y]+0.5)
|
||
|
img[y,x,:] = np.array(c)*255
|
||
|
#img[:,:,0] = np.minimum(np.abs(diff*255).astype(int), 255)[:Ny,:Nx]
|
||
|
out = img.reshape((Nx*Ny*3,)).astype(np.uint8)
|
||
|
#A += np.random.uniform(-0.01, 0.01, (N, N))
|
||
|
|
||
|
#print(len(out))
|
||
|
os.write(1, out.tobytes())
|
||
|
|
||
|
if solve and np.sum(np.abs(diff)) <= 0.001*N*N:
|
||
|
solve = False
|
||
|
elif not solve and np.sum(np.abs(diff)) >= 1*N*N:
|
||
|
solve = True
|
||
|
|
||
|
#os.write(2, b"frame")
|
||
|
#print(angles, diff)
|
||
|
|
||
|
time.sleep(0.01)
|