-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathcommon.py
105 lines (90 loc) · 2.37 KB
/
common.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
from numba import njit
from numpy import empty
from math import sqrt
from random import uniform
@njit
def vec2(x,y):
'''Convert two numbers to a 2D numpy array '''
v = empty((2,),dtype=float)
v[0] = x
v[1] = y
return v
@njit
def vec3(x,y,z):
'''Convert three numbers to a 2D numpy array '''
v = empty((3,),dtype=float)
v[0] = x
v[1] = y
v[2] = z
return v
@njit
def vec4(x,y,z,w):
'''Convert four numbers to a 2D numpy array '''
v = empty((4,),dtype=float)
v[0] = x
v[1] = y
v[2] = z
v[3] = w
return v
# @njit
# def V(x, inv_noise_scale, offset):
# ''' V returns a curl noise vector field'''
# p = empty(3)
# p[0:2] = x*inv_noise_scale+offset
# return scnoise_curl(p)[0:2]
@njit
def V(x, inv_noise_scale, offset, noise_grad):
''' V returns a curl noise vector field'''
p = empty(3)
p[0:2] = x*inv_noise_scale+offset
p[2] = 98.45
g = noise_grad(p)
gh = vec2(-g[1], g[0])
return gh
@njit
def Euler(x, dt, inv_noise_scale, offset, noise_grad):
'''Perform one Euler step of length dt from x along the curl noise vector field'''
return x + dt * V(x, inv_noise_scale, offset, noise_grad)
@njit
def RK4(x, dt, inv_noise_scale, offset, noise_grad):
'''Perform one Runge-Kutta 4th order step of length dt from x along the curl noise vector field'''
a = dt * V(x, inv_noise_scale, offset, noise_grad)
b = dt * V(x+a/2, inv_noise_scale, offset, noise_grad)
c = dt * V(x+b/2, inv_noise_scale, offset, noise_grad)
d = dt * V(x+c, inv_noise_scale, offset, noise_grad)
return x + (a+2*b+2*c+d)/6
@njit
def rand_vec2():
'''Convert two numbers to a 2D numpy array '''
v = empty((2,),dtype=float)
v[0] = uniform(0,1)
v[1] = uniform(0,1)
return v
## UNUSED
@njit
def smoothstep(a,b,x):
t = x-a/b-a
if t < 0: return 0
if 1 < t: return 1
return 3*t**2-2*t**3
@njit
def Rand(x, dt, inv_noise_scale, offset):
'''Perform a step of length dt from x in a random direction.'''
return x + dt * rand_vec2()
@njit
def NoiseV(x, dt, inv_noise_scale, offset, noise_grad):
'''Perform a step of length dt from x in a noise vector direction.'''
p = vec3(0,0,offset)
p[0:2] = x*inv_noise_scale
g = noise_grad(p)
return x + dt * g[0:2]
def L2_star_disc(P):
n = P.shape[0]
l2sd = (1/3)**2
l2sd -= (2/n) * sum([ 0.25 * (1 - x[0]**2) * (1 - x[1]**2) for x in P])
s = 0
for xi in P:
for xk in P:
s += (1-max(xi[0], xk[0]))*(1-max(xi[1], xk[1]))
l2sd += s/(n**2)
return sqrt(l2sd)