-
Notifications
You must be signed in to change notification settings - Fork 1
/
models.py
64 lines (59 loc) · 2.45 KB
/
models.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
import numpy as np
import random
import math
import scipy.stats
import scipy.constants
import loas
def get_Q_sfc(normal):
"""
Returns the quaternion to pass from the reference frame to
the local surface frame where the x axis is normal to the surface
"""
dir_rot = normal + loas.utils.tov(1,0,0)
if np.linalg.norm(dir_rot) < 1e-6:
dir_rot = loas.utils.tov(0,1,0)
return loas.utils.Quaternion(0, *dir_rot)
def maxwell(epsilon):
# See Sharipov, Rarefied gas dynamics, 4
def model(part_speed_i, normal, sat_temp, part_mass):
Q_sfc = get_Q_sfc(normal) #quaternion de passage sur la surface
if random.random() < epsilon:
# specular reflexion
normal_rel_speed = (np.transpose(normal) @ part_speed_i)[0,0]*normal
part_speed_r = part_speed_i - 2*normal_rel_speed
else:
part_speed_r_norm = scipy.stats.maxwell.rvs(scale = math.sqrt(scipy.constants.k*sat_temp/part_mass))
theta = math.asin(2*random.random()-1) #angle par rapport à la normale à la surface (donc le vecteur (1,0,0)) dans le repère de la sfc
phi = 2*math.pi*random.random() #angle dans le plan (yOz)
part_speed_r = Q_sfc.V2R(
part_speed_r_norm*
loas.utils.tov(
math.cos(theta),
math.sin(theta)*math.cos(phi),
math.sin(theta)*math.sin(phi)
)
)
return part_speed_r
return model
def schamberg(theta_, theta_i):
def model(part_speed_i, normal, sat_temp, part_mass):
# Schamberg model
Q_sfc = get_Q_sfc(normal)
part_speed_r_norm = scipy.stats.maxwell.rvs(scale = math.sqrt(scipy.constants.k*sat_temp/part_mass))
theta = 2*theta_0/math.pi*math.asin(2*random.random()-1) + theta_i #angle par rapport à la normale à la surface (donc le vecteur (1,0,0)) dans le repère de la sfc
phi = 2*math.pi*random.random() #angle dans le plan (yOz)
part_speed_r = Q_sfc.V2R(
part_speed_r_norm*
loas.utils.tov(
math.cos(theta),
math.sin(theta)*math.cos(phi),
math.sin(theta)*math.sin(phi)
)
)
return part_speed_r
return model
def schamberg_compose(schambergs, coefs):
def model(*args, **kwargs):
used_schamberg = random.choice(schambergs, coefs)
return used_schamberg(*args,**kwargs)
return model