clean slate

This commit is contained in:
zomseffen 2020-12-12 10:32:27 +01:00
parent 6577ac6558
commit afc0fff4fa

View file

@ -42,6 +42,7 @@ def value_to_color(v, min_value, max_value):
class Client: class Client:
def __init__(self, test=False, pos=[0, 0, 0]): def __init__(self, test=False, pos=[0, 0, 0]):
self.state = 0
with open('./config.json', 'r') as f: with open('./config.json', 'r') as f:
self.config = json.load(f) self.config = json.load(f)
glutInit(sys.argv) glutInit(sys.argv)
@ -146,6 +147,7 @@ class Client:
self.n_a_eq = np.zeros(self.n_a.shape) self.n_a_eq = np.zeros(self.n_a.shape)
self.n = np.zeros(self.field) self.n = np.zeros(self.field)
self.n[:, :, :] += 1.0 self.n[:, :, :] += 1.0
self.gravity_applies = np.zeros(self.field)
# self.n /= np.sum(self.n) # self.n /= np.sum(self.n)
self.n_a[0] = np.array(self.n) self.n_a[0] = np.array(self.n)
self.u = np.zeros(self.field + (self.e_a.shape[1],)) self.u = np.zeros(self.field + (self.e_a.shape[1],))
@ -214,134 +216,35 @@ class Client:
min_value = 0 min_value = 0
max_value_n = np.max(self.n) max_value_n = np.max(self.n)
# max_value = 1.0 # max_value_n = 1.0
vel = np.sqrt(np.sum(np.square(self.u), axis=3)) vel = np.sqrt(np.sum(np.square(self.u), axis=3)) *self.n
max_value_vel = np.max(vel) max_value_vel = np.max(vel)
# max_value_vel = np.sqrt(3) # max_value_vel = np.sqrt(3)
print('round') print('round')
print(max_value_n) print('sum n: %f' % np.sum(self.n))
print('max n: %f' % np.max(self.n))
print('min n: %f' % np.min(self.n))
print('sum vel: %f' % np.sum(vel))
print('max vel: %f' % np.max(vel))
print('min vel: %f' % np.min(vel))
for x_pos in range(0, 100): for x_pos in range(0, 100):
for y_pos in range(0, 100): for y_pos in range(0, 100):
for z_pos in range(0, 1): for z_pos in range(0, 1):
if self.state == 2:
r, g, b = value_to_color(int(self.gravity_applies[x_pos, y_pos, z_pos]), 0, 1)
if self.state == 1:
r, g, b = value_to_color(vel[x_pos, y_pos, z_pos], min_value, max_value_vel)
if self.state == 0:
r, g, b = value_to_color(self.n[x_pos, y_pos, z_pos], min_value, max_value_n) r, g, b = value_to_color(self.n[x_pos, y_pos, z_pos], min_value, max_value_n)
# r, g, b = value_to_color(vel[x_pos, y_pos, z_pos], min_value, max_value_vel)
self.world_provider.world.set_color(x_pos, y_pos, z_pos, r, g, b) self.world_provider.world.set_color(x_pos, y_pos, z_pos, r, g, b)
self.world_provider.world.set_color(int(round(self.test_pixel[0])), self.world_provider.world.set_color(int(round(self.test_pixel[0])),
int(round(self.test_pixel[1])), int(round(self.test_pixel[1])),
int(round(self.test_pixel[2])), 1.0, 1.0, 1.0) int(round(self.test_pixel[2])), 1.0, 1.0, 1.0)
# self.u *= 0.95
old_n_sum = np.sum(self.n)
for a in range(len(self.e_a)):
e_au = np.sum(self.e_a[a] * self.u, axis=3)
uu = np.sum(self.u * self.u, axis=3)
self.n_a_eq[a] = self.w_a[a] * self.n * (1.0 + 3.0 * e_au + 4.5 * np.square(e_au) - 1.5 * uu)
print(np.max(self.n_a_eq[0]))
if not self.compressible:
excess = (self.n_a_eq[0] > self.max_n) * (self.n_a_eq[0] - self.max_n)
dir_sum = np.sum(self.n_a_eq[1:], axis=0)
self.n_a_eq[1:] += excess * 1/8 #(self.n_a_eq[1:] / dir_sum)
self.n_a_eq[0] -= excess
n_a_t_1 = np.zeros((len(self.e_a),) + self.field)
for a in range(len(self.e_a)):
temp = ((-1.0 / self.relaxation_time) * (self.n_a[a] - self.n_a_eq[a]) + self.n_a[a])
n_a_t_1[a, max(self.e_a[a][0], 0): min(self.field[0], self.field[0] + self.e_a[a][0]),
max(self.e_a[a][1], 0): min(self.field[1], self.field[1] + self.e_a[a][1]),
max(self.e_a[a][2], 0): min(self.field[2], self.field[2] + self.e_a[a][2])] += temp[
max(-self.e_a[a][0], 0): min(
self.field[0],
self.field[0] -
self.e_a[a][0]),
max(-self.e_a[a][1], 0): min(
self.field[1],
self.field[1] -
self.e_a[a][1]),
max(-self.e_a[a][2], 0): min(
self.field[2],
self.field[2] -
self.e_a[a][2])]
for index in range(len(self.e_a[a])):
if self.e_a[a][index] != 0:
e_a_clipped = -np.array(self.e_a[a])
# e_a_clipped[index] = 0
# e_a_clipped[index] = -self.e_a[a][index]
clipped_dir = np.zeros(3)
clipped_dir[index] = self.e_a[a][index]
e_a_index = np.where(np.all(self.e_a == e_a_clipped, axis=1))[0][0]
# e_a_index = 0
if index == 0:
if self.e_a[a][index] > 0:
slice_index = -1
else:
slice_index = 0
n_a_t_1[e_a_index, slice_index, :, :] += temp[slice_index, :, :]
temp[slice_index, :, :] *= 0
if index == 1:
if self.e_a[a][index] > 0:
slice_index = -1
else:
slice_index = 0
n_a_t_1[e_a_index, :, slice_index, :] += temp[:, slice_index, :]
temp[:, slice_index, :] *= 0
if index == 2:
if self.e_a[a][index] > 0:
slice_index = -1
else:
slice_index = 0
n_a_t_1[e_a_index, :, :, slice_index] += temp[:, :, slice_index]
temp[:, :, slice_index] *= 0
self.n_a = n_a_t_1
if np.min(self.n_a) < 0:
test = 1
self.n_a = np.maximum(0, self.n_a)
self.n_a = old_n_sum * self.n_a / np.sum(self.n_a)
self.n = np.sum(self.n_a, axis=0, keepdims=False)
# self.n = np.sum(np.abs(n_a_t_1), axis=0, keepdims=False)
# stabilise the number because of rounding errors
self.n = old_n_sum * self.n / np.sum(self.n)
self.u *= 0
for a in range(len(self.e_a)):
self.u[:, :, :, 0] += self.n_a[a] * self.e_a[a][0]
self.u[:, :, :, 1] += self.n_a[a] * self.e_a[a][1]
self.u[:, :, :, 2] += self.n_a[a] * self.e_a[a][2]
self.u[:, :, :, 0] /= self.n
self.u[:, :, :, 1] /= self.n
self.u[:, :, :, 2] /= self.n
self.u[self.n == 0] = 0
length = np.sqrt(np.sum(np.square(self.u), axis=3, keepdims=True))
# gravity
gravity_applies = self.n < self.w_a[0]
gravity_applies[:, :-1, :] = gravity_applies[:, 1:, :]
gravity_applies[:, -1, :] = False
self.u[gravity_applies, 1] -= 0.01
new_lengths = np.sqrt(np.sum(np.square(self.u), axis=3, keepdims=True))
self.u = self.u / new_lengths * length
zero_length = (new_lengths == 0)
self.u[zero_length[:, :, :, 0], :] = 0
u = self.u[int(self.test_pixel[0]), int(self.test_pixel[1]), int(self.test_pixel[2])]
self.test_pixel[0] = max(0, min(self.field[0] - 1, self.test_pixel[0] + u[0]))
self.test_pixel[1] = max(0, min(self.field[1] - 1, self.test_pixel[1] + u[1]))
self.test_pixel[2] = max(0, min(self.field[2] - 1, self.test_pixel[2] + u[2]))
# print(1.0 / (time.time() - self.time)) # print(1.0 / (time.time() - self.time))
self.time = time.time() self.time = time.time()
glutPostRedisplay() glutPostRedisplay()
@ -378,6 +281,9 @@ class Client:
if key == b'e': if key == b'e':
self.opening += 0.25 self.opening += 0.25
if key == b'+':
self.state = (self.state +1) % 3
if key == b'r': if key == b'r':
print(self.cx, self.cy, self.opening) print(self.cx, self.cy, self.opening)
# glutPostRedisplay() # glutPostRedisplay()