Shamrock 2025.10.0
Astrophysical Code
Loading...
Searching...
No Matches
kernel_collision.py
1import numpy as np
2
3# from numba import njit
4
5
6# @njit
7def kconst(K0, u, v):
8 """
9 Function to compute constant kernel
10
11 Parameters
12 ----------
13 K0 : scalar, type -> float
14 constant value of the kernel function (used to adapt to code unit)
15 u : scalar, type -> float
16 mass variable (colliding grain of mass u)
17 v : scalar, type -> float
18 mass variable (colliding grain of mass v)
19
20
21 Returns
22 -------
23 res : scalar, type -> float
24 evaluate constant kernel
25
26 """
27
28 res = K0
29
30 return res
31
32
33# @njit
34def kadd(K0, u, v):
35 """
36 Function to compute additive kernel
37
38 Parameters
39 ----------
40 K0 : scalar, type -> float
41 constant value of the kernel function (used to adapt to code unit)
42 u : scalar, type -> float
43 mass variable (colliding grain of mass u)
44 v : scalar, type -> float
45 mass variable (colliding grain of mass v)
46
47
48 Returns
49 -------
50 res : scalar, type -> float
51 evaluate additive kernel at u and v
52
53 """
54
55 res = K0 * (u + v)
56 return res
57
58
59# @njit
60def kdv(K0, u, v):
61 """
62 Function to compute the cross-section term in the ballistic kernel K = sigma * dv
63
64 Parameters
65 ----------
66 K0 : scalar, type -> float
67 constant value of the kernel function (used to adapt to code unit)
68 u : scalar, type -> float
69 mass variable (colliding grain of mass u)
70 v : scalar, type -> float
71 mass variable (colliding grain of mass v)
72
73
74 Returns
75 -------
76 res : scalar, type -> float
77 evaluate cross-section term of the balistic kernel at u and v
78
79 """
80
81 # cross-section is calculated from mass variables
82 res = K0 * (u ** (2.0 / 3.0) + 2.0 * u ** (1.0 / 3.0) * v ** (1.0 / 3.0) + v ** (2.0 / 3.0))
83
84 return res
85
86
87# @njit
88def k_Br(K0, u, v):
89 """
90 Function to compute collision kernel from Brownian motion, K = sigma * dv
91
92 Parameters
93 ----------
94 K0 : scalar, type -> float
95 constant value of the kernel function (used to adapt to code unit)
96 u : scalar, type -> float
97 mass variable (colliding grain of mass u)
98 v : scalar, type -> float
99 mass variable (colliding grain of mass v)
100
101
102 Returns
103 -------
104 res : scalar, type -> float
105 Brownian motion collision kernel
106
107 """
108
109 res = (
110 K0
111 * (u ** (2.0 / 3.0) + 2.0 * u ** (1.0 / 3.0) * v ** (1.0 / 3.0) + v ** (2.0 / 3.0))
112 * np.sqrt(1.0 / u + 1.0 / v)
113 )
114
115 return res
116
117
118def func_kernel(kernel, K0, u, v):
119 """
120 Function to compute kernels at u and v
121
122 Parameters
123 ----------
124 kernel : scalar, type -> integer
125 select the collisional kernel function
126 K0 : scalar, type -> float
127 constant value of the kernel function (used to adapt to code unit)
128 u : scalar, type -> float
129 mass variable (colliding grain of mass u)
130 v : scalar, type -> float
131 mass variable (colliding grain of mass v)
132
133
134 Returns
135 -------
136 res : scalar, type -> float
137 evaluate kernel at u and v
138
139 """
140
141 match kernel:
142 case 0:
143 res = kconst(K0, u, v)
144 case 1:
145 res = kadd(K0, u, v)
146 case 2:
147 res = k_Br(K0, u, v)
148 case 3:
149 res = kdv(K0, u, v)
150 case _:
151 return "Need to choose a kernel in the list."
152
153 return res
154
155
156# @njit
157def func_kernel_numba(kernel, K0, u, v):
158 """
159 Function to compute kernels at u and v
160
161 Parameters
162 ----------
163 kernel : scalar, type -> integer
164 select the collisional kernel function
165 K0 : scalar, type -> float
166 constant value of the kernel function (used to adapt to code unit)
167 u : scalar, type -> float
168 mass variable (colliding grain of mass u)
169 v : scalar, type -> float
170 mass variable (colliding grain of mass v)
171
172
173 Returns
174 -------
175 res : scalar, type -> float
176 evaluate kernel at u and v
177
178 """
179
180 if kernel == 0:
181 res = kconst(K0, u, v)
182 elif kernel == 1:
183 res = kadd(K0, u, v)
184 elif kernel == 2:
185 res = k_Br(K0, u, v)
186 elif kernel == 3:
187 res = kdv(K0, u, v)
188 else:
189 # Return special value for unsupported kernel; cannot return string in njit
190 res = -1.0
191 return res